You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
4559 lines
112 KiB
4559 lines
112 KiB
/*=========================================================================
|
|
|
|
Program: Visualization Toolkit
|
|
Module: $RCSfile: vtkDistributedDataFilter.cxx,v $
|
|
|
|
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
|
|
All rights reserved.
|
|
See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
|
|
|
|
This software is distributed WITHOUT ANY WARRANTY; without even
|
|
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
|
|
PURPOSE. See the above copyright notice for more information.
|
|
|
|
=========================================================================*/
|
|
/*----------------------------------------------------------------------------
|
|
Copyright (c) Sandia Corporation
|
|
See Copyright.txt or http://www.paraview.org/HTML/Copyright.html for details.
|
|
----------------------------------------------------------------------------*/
|
|
|
|
// .NAME vtkDistributedDataFilter
|
|
//
|
|
// .SECTION Description
|
|
//
|
|
// .SECTION See Also
|
|
|
|
#include "vtkToolkits.h"
|
|
#include "vtkDistributedDataFilter.h"
|
|
#include "vtkModelMetadata.h"
|
|
#include "vtkExtractCells.h"
|
|
#include "vtkMergeCells.h"
|
|
#include "vtkInformation.h"
|
|
#include "vtkInformationVector.h"
|
|
#include "vtkStreamingDemandDrivenPipeline.h"
|
|
#include "vtkObjectFactory.h"
|
|
#include "vtkPKdTree.h"
|
|
#include "vtkUnstructuredGrid.h"
|
|
#include "vtkExtractUserDefinedPiece.h"
|
|
#include "vtkCellData.h"
|
|
#include "vtkCellArray.h"
|
|
#include "vtkPointData.h"
|
|
#include "vtkIntArray.h"
|
|
#include "vtkCharArray.h"
|
|
#include "vtkFloatArray.h"
|
|
#include "vtkUnsignedCharArray.h"
|
|
#include "vtkMultiProcessController.h"
|
|
#include "vtkSocketController.h"
|
|
#include "vtkDataSetWriter.h"
|
|
#include "vtkDataSetReader.h"
|
|
#include "vtkBoxClipDataSet.h"
|
|
#include "vtkClipDataSet.h"
|
|
#include "vtkBox.h"
|
|
#include "vtkIdList.h"
|
|
#include "vtkPointLocator.h"
|
|
#include "vtkPlane.h"
|
|
|
|
#ifdef VTK_USE_MPI
|
|
#include "vtkMPIController.h"
|
|
#endif
|
|
|
|
vtkCxxRevisionMacro(vtkDistributedDataFilter, "$Revision: 1.28 $")
|
|
|
|
vtkStandardNewMacro(vtkDistributedDataFilter)
|
|
|
|
#define TEMP_ELEMENT_ID_NAME "___D3___GlobalCellIds"
|
|
#define TEMP_INSIDE_BOX_FLAG "___D3___WHERE"
|
|
#define TEMP_NODE_ID_NAME "___D3___GlobalNodeIds"
|
|
|
|
#include <vtkstd/set>
|
|
#include <vtkstd/map>
|
|
#include <vtkstd/algorithm>
|
|
|
|
class vtkDistributedDataFilterSTLCloak
|
|
{
|
|
public:
|
|
vtkstd::map<int, int> IntMap;
|
|
vtkstd::multimap<int, int> IntMultiMap;
|
|
};
|
|
|
|
vtkDistributedDataFilter::vtkDistributedDataFilter()
|
|
{
|
|
this->Kdtree = NULL;
|
|
|
|
this->Controller = NULL;
|
|
this->SetController(vtkMultiProcessController::GetGlobalController());
|
|
|
|
this->Target = NULL;
|
|
this->Source = NULL;
|
|
|
|
this->NumConvexSubRegions = 0;
|
|
this->ConvexSubRegionBounds = NULL;
|
|
|
|
this->GhostLevel = 0;
|
|
|
|
this->GlobalNodeIdArrayName = NULL;
|
|
this->GlobalElementIdArrayName = NULL;
|
|
|
|
this->RetainKdtree = 1;
|
|
this->IncludeAllIntersectingCells = 0;
|
|
this->ClipCells = 0;
|
|
|
|
this->Timing = 0;
|
|
|
|
this->UseMinimalMemory = 0;
|
|
}
|
|
|
|
vtkDistributedDataFilter::~vtkDistributedDataFilter()
|
|
{
|
|
if (this->Kdtree)
|
|
{
|
|
this->Kdtree->Delete();
|
|
this->Kdtree = NULL;
|
|
}
|
|
|
|
this->SetController(NULL);
|
|
|
|
if (this->Target)
|
|
{
|
|
delete [] this->Target;
|
|
this->Target= NULL;
|
|
}
|
|
|
|
if (this->Source)
|
|
{
|
|
delete [] this->Source;
|
|
this->Source= NULL;
|
|
}
|
|
|
|
if (this->ConvexSubRegionBounds)
|
|
{
|
|
delete [] this->ConvexSubRegionBounds;
|
|
this->ConvexSubRegionBounds = NULL;
|
|
}
|
|
|
|
if (this->GlobalNodeIdArrayName)
|
|
{
|
|
delete [] this->GlobalNodeIdArrayName;
|
|
}
|
|
|
|
if (this->GlobalElementIdArrayName)
|
|
{
|
|
delete [] this->GlobalElementIdArrayName;
|
|
}
|
|
}
|
|
//-------------------------------------------------------------------------
|
|
// Global element and node IDs:
|
|
// Either the user gives us the names of these arrays, or we find them
|
|
// in the input (using the Exodus reader names for them), or we created
|
|
// these arrays.
|
|
//-------------------------------------------------------------------------
|
|
|
|
const char *vtkDistributedDataFilter::GetGlobalElementIdArrayName(vtkDataSet *set)
|
|
{
|
|
//------------------------------------------------
|
|
// list common names for global element id arrays here
|
|
//
|
|
int nnames = 1;
|
|
const char *arrayNames[1] = {
|
|
"GlobalElementId" // vtkExodusReader name
|
|
};
|
|
//------------------------------------------------
|
|
|
|
// ParaView does this... we need to fix it.
|
|
if (this->GlobalElementIdArrayName && (!this->GlobalElementIdArrayName[0]))
|
|
{
|
|
delete [] this->GlobalElementIdArrayName;
|
|
this->GlobalElementIdArrayName = NULL;
|
|
}
|
|
|
|
const char *gidArrayName = NULL;
|
|
|
|
if (this->GlobalElementIdArrayName)
|
|
{
|
|
vtkDataArray *da = set->GetCellData()->GetArray(this->GlobalElementIdArrayName);
|
|
|
|
if (da)
|
|
{
|
|
gidArrayName = this->GlobalElementIdArrayName;
|
|
}
|
|
else
|
|
{
|
|
this->SetGlobalElementIdArrayName(NULL);
|
|
}
|
|
}
|
|
|
|
if (!gidArrayName)
|
|
{
|
|
// Maybe we can find a global element ID array
|
|
|
|
for (int nameId=0; nameId < nnames; nameId++)
|
|
{
|
|
vtkDataArray *da = set->GetCellData()->GetArray(arrayNames[nameId]);
|
|
|
|
if (da)
|
|
{
|
|
this->SetGlobalElementIdArrayName(arrayNames[nameId]);
|
|
gidArrayName = arrayNames[nameId];
|
|
break;
|
|
}
|
|
}
|
|
}
|
|
|
|
return gidArrayName;
|
|
}
|
|
int *vtkDistributedDataFilter::GetGlobalElementIds(vtkDataSet *set)
|
|
{
|
|
const char *geidName = this->GetGlobalElementIdArrayName(set);
|
|
|
|
if (!geidName)
|
|
{
|
|
return NULL;
|
|
}
|
|
|
|
vtkDataArray *da = set->GetCellData()->GetArray(geidName);
|
|
vtkIntArray *ia = vtkIntArray::SafeDownCast(da);
|
|
|
|
if (!ia)
|
|
{
|
|
return NULL;
|
|
}
|
|
|
|
return ia->GetPointer(0);
|
|
}
|
|
|
|
const char *vtkDistributedDataFilter::GetGlobalNodeIdArrayName(vtkDataSet *set)
|
|
{
|
|
//------------------------------------------------
|
|
// list common names for global node id arrays here
|
|
//
|
|
int nnames = 1;
|
|
const char *arrayNames[1] = {
|
|
"GlobalNodeId" // vtkExodusReader name
|
|
};
|
|
//------------------------------------------------
|
|
|
|
// ParaView does this... we need to fix it.
|
|
if (this->GlobalNodeIdArrayName && (!this->GlobalNodeIdArrayName[0]))
|
|
{
|
|
delete [] this->GlobalNodeIdArrayName;
|
|
this->GlobalNodeIdArrayName = NULL;
|
|
}
|
|
|
|
const char *gidArrayName = NULL;
|
|
|
|
if (this->GlobalNodeIdArrayName)
|
|
{
|
|
vtkDataArray *da = set->GetPointData()->GetArray(this->GlobalNodeIdArrayName);
|
|
|
|
if (da)
|
|
{
|
|
gidArrayName = this->GlobalNodeIdArrayName;
|
|
}
|
|
else
|
|
{
|
|
this->SetGlobalNodeIdArrayName(NULL);
|
|
}
|
|
}
|
|
|
|
if (!gidArrayName)
|
|
{
|
|
// Maybe we can find a global node ID array
|
|
|
|
for (int nameId=0; nameId < nnames; nameId++)
|
|
{
|
|
vtkDataArray *da = set->GetPointData()->GetArray(arrayNames[nameId]);
|
|
|
|
if (da)
|
|
{
|
|
this->SetGlobalNodeIdArrayName(arrayNames[nameId]);
|
|
gidArrayName = arrayNames[nameId];
|
|
break;
|
|
}
|
|
}
|
|
}
|
|
|
|
return gidArrayName;
|
|
}
|
|
int *vtkDistributedDataFilter::GetGlobalNodeIds(vtkDataSet *set)
|
|
{
|
|
const char *gnidName = this->GetGlobalNodeIdArrayName(set);
|
|
|
|
if (!gnidName)
|
|
{
|
|
return NULL;
|
|
}
|
|
|
|
vtkDataArray *da = set->GetPointData()->GetArray(gnidName);
|
|
vtkIntArray *ia = vtkIntArray::SafeDownCast(da);
|
|
|
|
if (!ia)
|
|
{
|
|
return NULL;
|
|
}
|
|
|
|
return ia->GetPointer(0);
|
|
}
|
|
|
|
//-------------------------------------------------------------------------
|
|
//-------------------------------------------------------------------------
|
|
vtkPKdTree *vtkDistributedDataFilter::GetKdtree()
|
|
{
|
|
if (this->Kdtree == NULL)
|
|
{
|
|
this->Kdtree = vtkPKdTree::New();
|
|
this->Kdtree->AssignRegionsContiguous();
|
|
this->Kdtree->SetTiming(this->GetTiming());
|
|
}
|
|
|
|
return this->Kdtree;
|
|
}
|
|
unsigned long vtkDistributedDataFilter::GetMTime()
|
|
{
|
|
unsigned long t1, t2;
|
|
|
|
t1 = this->Superclass::GetMTime();
|
|
if (this->Kdtree == NULL)
|
|
{
|
|
return t1;
|
|
}
|
|
t2 = this->Kdtree->GetMTime();
|
|
if (t1 > t2)
|
|
{
|
|
return t1;
|
|
}
|
|
return t2;
|
|
}
|
|
|
|
void vtkDistributedDataFilter::SetController(vtkMultiProcessController *c)
|
|
{
|
|
if (this->Kdtree)
|
|
{
|
|
this->Kdtree->SetController(c);
|
|
}
|
|
|
|
if ((c == NULL) || (c->GetNumberOfProcesses() == 0))
|
|
{
|
|
this->NumProcesses = 1;
|
|
this->MyId = 0;
|
|
}
|
|
|
|
if (this->Controller == c)
|
|
{
|
|
return;
|
|
}
|
|
|
|
this->Modified();
|
|
|
|
if (this->Controller != NULL)
|
|
{
|
|
this->Controller->UnRegister(this);
|
|
this->Controller = NULL;
|
|
}
|
|
|
|
if (c == NULL)
|
|
{
|
|
return;
|
|
}
|
|
|
|
this->Controller = c;
|
|
|
|
c->Register(this);
|
|
this->NumProcesses = c->GetNumberOfProcesses();
|
|
this->MyId = c->GetLocalProcessId();
|
|
}
|
|
|
|
void vtkDistributedDataFilter::SetBoundaryMode(int mode)
|
|
{
|
|
switch (mode)
|
|
{
|
|
case vtkDistributedDataFilter::ASSIGN_TO_ONE_REGION:
|
|
this->AssignBoundaryCellsToOneRegionOn();
|
|
break;
|
|
case vtkDistributedDataFilter::ASSIGN_TO_ALL_INTERSECTING_REGIONS:
|
|
this->AssignBoundaryCellsToAllIntersectingRegionsOn();
|
|
break;
|
|
case vtkDistributedDataFilter::SPLIT_BOUNDARY_CELLS:
|
|
this->DivideBoundaryCellsOn();
|
|
break;
|
|
}
|
|
}
|
|
|
|
int vtkDistributedDataFilter::GetBoundaryMode()
|
|
{
|
|
if (!this->IncludeAllIntersectingCells && !this->ClipCells)
|
|
{
|
|
return vtkDistributedDataFilter::ASSIGN_TO_ONE_REGION;
|
|
}
|
|
if (this->IncludeAllIntersectingCells && !this->ClipCells)
|
|
{
|
|
return vtkDistributedDataFilter::ASSIGN_TO_ALL_INTERSECTING_REGIONS;
|
|
}
|
|
if (this->IncludeAllIntersectingCells && this->ClipCells)
|
|
{
|
|
return vtkDistributedDataFilter::SPLIT_BOUNDARY_CELLS;
|
|
}
|
|
|
|
return -1;
|
|
}
|
|
|
|
void vtkDistributedDataFilter::AssignBoundaryCellsToOneRegionOn()
|
|
{
|
|
this->SetAssignBoundaryCellsToOneRegion(1);
|
|
}
|
|
void vtkDistributedDataFilter::AssignBoundaryCellsToOneRegionOff()
|
|
{
|
|
this->SetAssignBoundaryCellsToOneRegion(0);
|
|
}
|
|
void vtkDistributedDataFilter::SetAssignBoundaryCellsToOneRegion(int val)
|
|
{
|
|
if (val)
|
|
{
|
|
this->IncludeAllIntersectingCells = 0;
|
|
this->ClipCells = 0;
|
|
}
|
|
}
|
|
void
|
|
vtkDistributedDataFilter::AssignBoundaryCellsToAllIntersectingRegionsOn()
|
|
{
|
|
this->SetAssignBoundaryCellsToAllIntersectingRegions(1);
|
|
}
|
|
void
|
|
vtkDistributedDataFilter::AssignBoundaryCellsToAllIntersectingRegionsOff()
|
|
{
|
|
this->SetAssignBoundaryCellsToAllIntersectingRegions(0);
|
|
}
|
|
void
|
|
vtkDistributedDataFilter::SetAssignBoundaryCellsToAllIntersectingRegions(int val)
|
|
{
|
|
if (val)
|
|
{
|
|
this->IncludeAllIntersectingCells = 1;
|
|
this->ClipCells = 0;
|
|
}
|
|
}
|
|
void vtkDistributedDataFilter::DivideBoundaryCellsOn()
|
|
{
|
|
this->SetDivideBoundaryCells(1);
|
|
}
|
|
void vtkDistributedDataFilter::DivideBoundaryCellsOff()
|
|
{
|
|
this->SetDivideBoundaryCells(0);
|
|
}
|
|
void vtkDistributedDataFilter::SetDivideBoundaryCells(int val)
|
|
{
|
|
if (val)
|
|
{
|
|
this->IncludeAllIntersectingCells = 1;
|
|
this->ClipCells = 1;
|
|
}
|
|
}
|
|
//-------------------------------------------------------------------------
|
|
// Execute
|
|
//-------------------------------------------------------------------------
|
|
|
|
int vtkDistributedDataFilter::RequestUpdateExtent(
|
|
vtkInformation *vtkNotUsed(request),
|
|
vtkInformationVector **inputVector,
|
|
vtkInformationVector *outputVector)
|
|
{
|
|
// get the info objects
|
|
vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
|
|
vtkInformation *outInfo = outputVector->GetInformationObject(0);
|
|
|
|
int piece, numPieces, ghostLevels;
|
|
|
|
// We require preceding filters to refrain from creating ghost cells.
|
|
|
|
piece = outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER());
|
|
numPieces = outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES());
|
|
ghostLevels = 0;
|
|
|
|
inInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER(), piece);
|
|
inInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES(),
|
|
numPieces);
|
|
inInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS(),
|
|
ghostLevels);
|
|
inInfo->Set(vtkStreamingDemandDrivenPipeline::EXACT_EXTENT(), 1);
|
|
|
|
return 1;
|
|
}
|
|
|
|
int vtkDistributedDataFilter::RequestInformation(
|
|
vtkInformation *vtkNotUsed(request),
|
|
vtkInformationVector **inputVector,
|
|
vtkInformationVector *outputVector)
|
|
{
|
|
// get the info objects
|
|
vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
|
|
vtkInformation *outInfo = outputVector->GetInformationObject(0);
|
|
|
|
outInfo->Set(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT(),
|
|
inInfo->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT()),
|
|
6);
|
|
outInfo->Set(vtkStreamingDemandDrivenPipeline::EXTENT_TRANSLATOR(),
|
|
inInfo->Get(vtkStreamingDemandDrivenPipeline::EXTENT_TRANSLATOR()));
|
|
outInfo->Set(vtkStreamingDemandDrivenPipeline::MAXIMUM_NUMBER_OF_PIECES(), -1);
|
|
|
|
return 1;
|
|
}
|
|
|
|
int vtkDistributedDataFilter::RequestData(
|
|
vtkInformation *vtkNotUsed(request),
|
|
vtkInformationVector **inputVector,
|
|
vtkInformationVector *outputVector)
|
|
{
|
|
// get the info objects
|
|
vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
|
|
vtkInformation *outInfo = outputVector->GetInformationObject(0);
|
|
|
|
// get the input and ouptut
|
|
vtkDataSet *input = vtkDataSet::SafeDownCast(
|
|
inInfo->Get(vtkDataObject::DATA_OBJECT()));
|
|
vtkUnstructuredGrid *output = vtkUnstructuredGrid::SafeDownCast(
|
|
outInfo->Get(vtkDataObject::DATA_OBJECT()));
|
|
|
|
this->GhostLevel = outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS());
|
|
this->NextProgressStep = 0;
|
|
int progressSteps = 5 + this->GhostLevel;
|
|
if (this->ClipCells)
|
|
{
|
|
progressSteps++;
|
|
}
|
|
|
|
this->ProgressIncrement = 1.0 / (double)progressSteps;
|
|
|
|
this->UpdateProgress(this->NextProgressStep++ * this->ProgressIncrement);
|
|
this->SetProgressText("Begin data redistribution");
|
|
|
|
if (this->NumProcesses == 1)
|
|
{
|
|
this->SingleProcessExecute(input, output);
|
|
this->UpdateProgress(1.0);
|
|
return 1;
|
|
}
|
|
|
|
// This method requires an MPI controller.
|
|
|
|
int aok = 0;
|
|
|
|
#ifdef VTK_USE_MPI
|
|
if (vtkMPIController::SafeDownCast(this->Controller))
|
|
{
|
|
aok = 1;
|
|
}
|
|
#endif
|
|
|
|
if (!aok)
|
|
{
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter multiprocess requires MPI");
|
|
return 1;
|
|
}
|
|
|
|
// Stage (0) - If any processes have 0 cell input data sets, then
|
|
// spread the input data sets around (quickly) before formal
|
|
// redistribution.
|
|
|
|
vtkDataSet *splitInput = this->TestFixTooFewInputFiles(input);
|
|
|
|
if (splitInput == NULL)
|
|
{
|
|
return 1; // Fewer cells than processes - can't divide input
|
|
}
|
|
|
|
this->UpdateProgress(this->NextProgressStep++ * this->ProgressIncrement);
|
|
this->SetProgressText("Compute spatial partitioning");
|
|
|
|
if (this->ClipCells && vtkDistributedDataFilter::HasMetadata(splitInput))
|
|
{
|
|
// Clipping cells invalidates metadata that is cell based
|
|
// Here we should remove the metadata and display a warning
|
|
}
|
|
|
|
// Stage (1) - use vtkPKdTree to...
|
|
// Create a load balanced spatial decomposition in parallel.
|
|
// Create a table assigning regions to processes.
|
|
//
|
|
// Note k-d tree will only be re-built if input or parameters
|
|
// have changed on any of the processing nodes.
|
|
|
|
int fail = this->PartitionDataAndAssignToProcesses(splitInput);
|
|
|
|
if (fail)
|
|
{
|
|
if (splitInput != input)
|
|
{
|
|
splitInput->Delete();
|
|
}
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::Execute k-d tree failure");
|
|
return 1;
|
|
}
|
|
|
|
this->UpdateProgress(this->NextProgressStep++ * this->ProgressIncrement);
|
|
this->SetProgressText("Compute global data array bounds");
|
|
|
|
// Let the vtkPKdTree class compile global bounds for all
|
|
// data arrays. These can be accessed by D3 user by getting
|
|
// a handle to the vtkPKdTree object and querying it.
|
|
|
|
this->Kdtree->CreateGlobalDataArrayBounds();
|
|
|
|
this->UpdateProgress(this->NextProgressStep++ * this->ProgressIncrement);
|
|
this->SetProgressText("Redistribute data");
|
|
|
|
// Stage (2) - Redistribute data, so that each process gets a ugrid
|
|
// containing the cells in it's assigned spatial regions. (Note
|
|
// that a side effect of merging the grids received from different
|
|
// processes is that the final grid has no duplicate points.)
|
|
//
|
|
// This call will delete splitInput if it's not this->GetInput().
|
|
|
|
vtkUnstructuredGrid *redistributedInput = this->RedistributeDataSet(splitInput,
|
|
input);
|
|
|
|
if (redistributedInput == NULL)
|
|
{
|
|
this->Kdtree->Delete();
|
|
this->Kdtree = NULL;
|
|
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::Execute redistribute failure");
|
|
return 1;
|
|
}
|
|
|
|
this->UpdateProgress(this->NextProgressStep++ * this->ProgressIncrement);
|
|
|
|
// Stage (3) - Add ghost cells to my sub grid.
|
|
|
|
vtkUnstructuredGrid *expandedGrid = redistributedInput;
|
|
|
|
if (this->GhostLevel > 0)
|
|
{
|
|
// Create global nodes IDs if we don't have them
|
|
|
|
if (this->GetGlobalNodeIdArrayName(redistributedInput) == NULL)
|
|
{
|
|
this->SetProgressText("Assign global point IDs");
|
|
int rc = this->AssignGlobalNodeIds(redistributedInput);
|
|
if (rc)
|
|
{
|
|
redistributedInput->Delete();
|
|
this->Kdtree->Delete();
|
|
this->Kdtree = NULL;
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::Execute global node id creation");
|
|
return 1;
|
|
}
|
|
}
|
|
|
|
// redistributedInput will be deleted by AcquireGhostCells
|
|
|
|
this->SetProgressText("Exchange ghost cells");
|
|
expandedGrid = this->AcquireGhostCells(redistributedInput);
|
|
}
|
|
|
|
// Stage (4) - Clip cells to the spatial region boundaries
|
|
|
|
if (this->ClipCells)
|
|
{
|
|
this->SetProgressText("Clip boundary cells");
|
|
this->ClipGridCells(expandedGrid);
|
|
this->UpdateProgress(this->NextProgressStep++ * this->ProgressIncrement);
|
|
}
|
|
|
|
// remove temporary arrays we created
|
|
|
|
this->SetProgressText("Clean up and finish");
|
|
|
|
vtkDataArray *da = expandedGrid->GetCellData()->GetArray(TEMP_ELEMENT_ID_NAME);
|
|
|
|
if (da)
|
|
{
|
|
expandedGrid->GetCellData()->RemoveArray(TEMP_ELEMENT_ID_NAME);
|
|
}
|
|
|
|
da = expandedGrid->GetPointData()->GetArray(TEMP_NODE_ID_NAME);
|
|
|
|
if (da)
|
|
{
|
|
expandedGrid->GetCellData()->RemoveArray(TEMP_NODE_ID_NAME);
|
|
}
|
|
|
|
output->ShallowCopy(expandedGrid);
|
|
|
|
expandedGrid->Delete();
|
|
|
|
if (!this->RetainKdtree)
|
|
{
|
|
this->Kdtree->Delete();
|
|
this->Kdtree = NULL;
|
|
}
|
|
else
|
|
{
|
|
this->Kdtree->SetDataSet(NULL);
|
|
}
|
|
|
|
this->UpdateProgress(1);
|
|
|
|
return 1;
|
|
}
|
|
vtkUnstructuredGrid *vtkDistributedDataFilter::RedistributeDataSet(
|
|
vtkDataSet *set, vtkDataSet *input)
|
|
{
|
|
// Create global cell ids before redistributing data. These
|
|
// will be necessary if we need ghost cells later on.
|
|
|
|
vtkDataSet *inputPlus = set;
|
|
|
|
if ((this->GhostLevel > 0) && (this->GetGlobalElementIdArrayName(set) == NULL))
|
|
{
|
|
if (set == input)
|
|
{
|
|
inputPlus = set->NewInstance();
|
|
inputPlus->ShallowCopy(set);
|
|
}
|
|
|
|
this->AssignGlobalElementIds(inputPlus);
|
|
}
|
|
|
|
// next call deletes inputPlus at the earliest opportunity
|
|
|
|
vtkUnstructuredGrid *finalGrid = this->MPIRedistribute(inputPlus, input);
|
|
|
|
return finalGrid;
|
|
}
|
|
|
|
int vtkDistributedDataFilter::PartitionDataAndAssignToProcesses(vtkDataSet *set)
|
|
{
|
|
if (this->Kdtree == NULL)
|
|
{
|
|
this->Kdtree = vtkPKdTree::New();
|
|
this->Kdtree->AssignRegionsContiguous();
|
|
this->Kdtree->SetTiming(this->GetTiming());
|
|
}
|
|
|
|
this->Kdtree->SetController(this->Controller);
|
|
this->Kdtree->SetNumberOfRegionsOrMore(this->NumProcesses);
|
|
this->Kdtree->SetMinCells(2);
|
|
|
|
this->Kdtree->SetDataSet(set);
|
|
|
|
// BuildLocator is smart enough to rebuild the k-d tree only if
|
|
// the input geometry has changed, or the k-d tree build parameters
|
|
// have changed. It will reassign regions if the region assignment
|
|
// scheme has changed.
|
|
|
|
this->Kdtree->BuildLocator();
|
|
|
|
int nregions = this->Kdtree->GetNumberOfRegions();
|
|
|
|
if (nregions < this->NumProcesses)
|
|
{
|
|
if (nregions == 0)
|
|
{
|
|
vtkErrorMacro("Unable to build k-d tree structure");
|
|
}
|
|
else
|
|
{
|
|
vtkErrorMacro("K-d tree must have at least one region per process");
|
|
}
|
|
this->Kdtree->Delete();
|
|
this->Kdtree = NULL;
|
|
return 1;
|
|
}
|
|
|
|
return 0;
|
|
}
|
|
int vtkDistributedDataFilter::ClipGridCells(vtkUnstructuredGrid *grid)
|
|
{
|
|
if (grid->GetNumberOfCells() == 0)
|
|
{
|
|
return 0;
|
|
}
|
|
|
|
// Global point IDs are meaningless after
|
|
// clipping, since this tetrahedralizes the whole data set.
|
|
// We remove that array.
|
|
|
|
const char *nodeIds = this->GetGlobalNodeIdArrayName(grid);
|
|
|
|
if (nodeIds)
|
|
{
|
|
grid->GetPointData()->RemoveArray(nodeIds);
|
|
this->GlobalNodeIdArrayName = NULL;
|
|
}
|
|
|
|
this->ClipCellsToSpatialRegion(grid);
|
|
|
|
return 0;
|
|
}
|
|
vtkUnstructuredGrid *
|
|
vtkDistributedDataFilter::AcquireGhostCells(vtkUnstructuredGrid *grid)
|
|
{
|
|
if (this->GhostLevel < 1)
|
|
{
|
|
return grid;
|
|
}
|
|
|
|
// Create a search structure mapping global point IDs to local point IDs
|
|
|
|
vtkIdType numPoints = grid->GetNumberOfPoints();
|
|
|
|
int *gnids = NULL;
|
|
|
|
if (numPoints > 0)
|
|
{
|
|
gnids = this->GetGlobalNodeIds(grid);
|
|
|
|
if (!gnids)
|
|
{
|
|
vtkWarningMacro(<< "Can't create ghost cells without global node IDs");
|
|
return grid;
|
|
}
|
|
}
|
|
|
|
vtkDistributedDataFilterSTLCloak *globalToLocalMap
|
|
= new vtkDistributedDataFilterSTLCloak;
|
|
|
|
for (int localPtId = 0; localPtId < numPoints; localPtId++)
|
|
{
|
|
const int id = gnids[localPtId];
|
|
globalToLocalMap->IntMap.insert(vtkstd::pair<const int, int>(id, localPtId));
|
|
}
|
|
|
|
vtkUnstructuredGrid *expandedGrid= NULL;
|
|
|
|
if (this->IncludeAllIntersectingCells)
|
|
{
|
|
expandedGrid =
|
|
this->AddGhostCellsDuplicateCellAssignment(grid, globalToLocalMap);
|
|
}
|
|
else
|
|
{
|
|
expandedGrid =
|
|
this->AddGhostCellsUniqueCellAssignment(grid, globalToLocalMap);
|
|
}
|
|
|
|
return expandedGrid;
|
|
}
|
|
void vtkDistributedDataFilter::SingleProcessExecute(vtkDataSet *input,
|
|
vtkUnstructuredGrid *output)
|
|
{
|
|
vtkDebugMacro(<< "vtkDistributedDataFilter::SingleProcessExecute()");
|
|
|
|
// we run the input through vtkMergeCells which will remove
|
|
// duplicate points
|
|
|
|
vtkDataSet* tmp = input->NewInstance();
|
|
tmp->ShallowCopy(input);
|
|
|
|
float tolerance = 0.0;
|
|
|
|
if (this->RetainKdtree)
|
|
{
|
|
if (this->Kdtree == NULL)
|
|
{
|
|
this->Kdtree = vtkPKdTree::New();
|
|
this->Kdtree->SetTiming(this->GetTiming());
|
|
}
|
|
|
|
this->Kdtree->SetDataSet(tmp);
|
|
this->Kdtree->BuildLocator();
|
|
tolerance = (float)this->Kdtree->GetFudgeFactor();
|
|
this->Kdtree->CreateGlobalDataArrayBounds();
|
|
}
|
|
else if (this->Kdtree)
|
|
{
|
|
this->Kdtree->Delete();
|
|
this->Kdtree = NULL;
|
|
}
|
|
|
|
vtkUnstructuredGrid *clean =
|
|
vtkDistributedDataFilter::MergeGrids(&tmp, 1, DeleteYes,
|
|
this->GetGlobalNodeIdArrayName(input), tolerance, NULL);
|
|
|
|
output->ShallowCopy(clean);
|
|
clean->Delete();
|
|
|
|
if (this->GhostLevel > 0)
|
|
{
|
|
// Add the vtkGhostLevels arrays. We have the whole
|
|
// data set, so all cells are level 0.
|
|
|
|
vtkDistributedDataFilter::AddConstantUnsignedCharPointArray(
|
|
output, "vtkGhostLevels", 0);
|
|
vtkDistributedDataFilter::AddConstantUnsignedCharCellArray(
|
|
output, "vtkGhostLevels", 0);
|
|
}
|
|
}
|
|
void vtkDistributedDataFilter::ComputeMyRegionBounds()
|
|
{
|
|
vtkIntArray *myRegions = vtkIntArray::New();
|
|
|
|
this->Kdtree->GetRegionAssignmentList(this->MyId, myRegions);
|
|
|
|
if (myRegions->GetNumberOfTuples() > 0)
|
|
{
|
|
this->NumConvexSubRegions =
|
|
this->Kdtree->MinimalNumberOfConvexSubRegions(
|
|
myRegions, &this->ConvexSubRegionBounds);
|
|
}
|
|
else
|
|
{
|
|
this->NumConvexSubRegions = 0;
|
|
if (this->ConvexSubRegionBounds)
|
|
{
|
|
delete [] this->ConvexSubRegionBounds;
|
|
this->ConvexSubRegionBounds = NULL;
|
|
}
|
|
}
|
|
|
|
myRegions->Delete();
|
|
}
|
|
int vtkDistributedDataFilter::CheckFieldArrayTypes(vtkDataSet *set)
|
|
{
|
|
int i;
|
|
|
|
// problem - vtkIdType arrays are written out as int arrays
|
|
// when marshalled with vtkDataWriter. This is a problem
|
|
// when receive the array and try to merge it with our own,
|
|
// which is a vtkIdType
|
|
|
|
vtkPointData *pd = set->GetPointData();
|
|
vtkCellData *cd = set->GetCellData();
|
|
|
|
int npointArrays = pd->GetNumberOfArrays();
|
|
|
|
for (i=0; i<npointArrays; i++)
|
|
{
|
|
int arrayType = pd->GetArray(i)->GetDataType();
|
|
|
|
if (arrayType == VTK_ID_TYPE)
|
|
{
|
|
return 1;
|
|
}
|
|
}
|
|
|
|
int ncellArrays = cd->GetNumberOfArrays();
|
|
|
|
for (i=0; i<ncellArrays; i++)
|
|
{
|
|
int arrayType = cd->GetArray(i)->GetDataType();
|
|
|
|
if (arrayType == VTK_ID_TYPE)
|
|
{
|
|
return 1;
|
|
}
|
|
}
|
|
|
|
return 0;
|
|
}
|
|
//-------------------------------------------------------------------------
|
|
// Quickly spread input data around if there are more processes than
|
|
// input data sets.
|
|
//-------------------------------------------------------------------------
|
|
extern "C"
|
|
{
|
|
int vtkDistributedDataFilterSortSize(const void *s1, const void *s2)
|
|
{
|
|
struct _procInfo{ int had; int procId; int has; } *a, *b;
|
|
|
|
a = (struct _procInfo *)s1;
|
|
b = (struct _procInfo *)s2;
|
|
|
|
if (a->has < b->has)
|
|
{
|
|
return 1;
|
|
}
|
|
else if (a->has == b->has)
|
|
{
|
|
return 0;
|
|
}
|
|
else
|
|
{
|
|
return -1;
|
|
}
|
|
}
|
|
}
|
|
vtkDataSet *vtkDistributedDataFilter::TestFixTooFewInputFiles(vtkDataSet *input)
|
|
{
|
|
int i, proc;
|
|
int me = this->MyId;
|
|
int nprocs = this->NumProcesses;
|
|
|
|
int numMyCells = input->GetNumberOfCells();
|
|
|
|
// Find out how many input cells each process has.
|
|
|
|
vtkIntArray *inputSize = this->ExchangeCounts(numMyCells, 0x0001);
|
|
int *sizes = inputSize->GetPointer(0);
|
|
|
|
int *nodeType = new int [nprocs];
|
|
const int Producer = 1;
|
|
const int Consumer = 2;
|
|
int numConsumers = 0;
|
|
int numTotalCells = 0;
|
|
|
|
for (proc = 0; proc < nprocs ; proc++)
|
|
{
|
|
numTotalCells += sizes[proc];
|
|
if (sizes[proc] == 0)
|
|
{
|
|
numConsumers++;
|
|
nodeType[proc] = Consumer;
|
|
}
|
|
else
|
|
{
|
|
nodeType[proc] = Producer;
|
|
}
|
|
}
|
|
|
|
if (numConsumers == 0)
|
|
{
|
|
// Nothing to do. Every process has input data.
|
|
|
|
delete [] nodeType;
|
|
inputSize->Delete();
|
|
return input;
|
|
}
|
|
|
|
if (numTotalCells < nprocs)
|
|
{
|
|
vtkErrorMacro(<< "D3 - fewer cells than processes");
|
|
delete [] nodeType;
|
|
inputSize->Delete();
|
|
return NULL;
|
|
}
|
|
|
|
int cellsPerNode = numTotalCells / nprocs;
|
|
|
|
vtkIdList **sendCells = new vtkIdList * [ nprocs ];
|
|
memset(sendCells, 0, sizeof(vtkIdList *) * nprocs);
|
|
|
|
if (numConsumers == nprocs - 1)
|
|
{
|
|
// Simple and common case.
|
|
// Only one process has data and divides it among the rest.
|
|
|
|
inputSize->Delete();
|
|
|
|
if (nodeType[me] == Producer)
|
|
{
|
|
int sizeLast = numTotalCells - ((nprocs-1) * cellsPerNode);
|
|
vtkIdType cellId = 0;
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
int ncells = ((proc == nprocs - 1) ? sizeLast :cellsPerNode);
|
|
|
|
sendCells[proc] = vtkIdList::New();
|
|
sendCells[proc]->SetNumberOfIds(ncells);
|
|
|
|
for (i=0; i<ncells; i++)
|
|
{
|
|
sendCells[proc]->SetId(i, cellId++);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
else
|
|
{
|
|
// The processes with data send it to processes without data.
|
|
// This is not the most balanced decomposition, and it is not the
|
|
// fastest. It is somewhere inbetween.
|
|
|
|
int minCells = (int)(.8 * cellsPerNode);
|
|
|
|
struct _procInfo{ int had; int procId; int has; };
|
|
|
|
struct _procInfo *procInfo = new struct _procInfo [nprocs];
|
|
|
|
for (proc = 0; proc < nprocs ; proc++)
|
|
{
|
|
procInfo[proc].had = inputSize->GetValue(proc);
|
|
procInfo[proc].procId = proc;
|
|
procInfo[proc].has = inputSize->GetValue(proc);
|
|
}
|
|
|
|
inputSize->Delete();
|
|
|
|
qsort(procInfo, nprocs, sizeof(struct _procInfo),
|
|
vtkDistributedDataFilterSortSize);
|
|
|
|
struct _procInfo *nextProducer = procInfo;
|
|
struct _procInfo *nextConsumer = procInfo + (nprocs - 1);
|
|
|
|
int numTransferCells = 0;
|
|
|
|
int sanityCheck=0;
|
|
int nprocsSquared = nprocs * nprocs;
|
|
|
|
while (sanityCheck++ < nprocsSquared)
|
|
{
|
|
int c = nextConsumer->procId;
|
|
|
|
if (nodeType[c] == Producer)
|
|
{
|
|
break;
|
|
}
|
|
|
|
int cGetMin = minCells - nextConsumer->has;
|
|
|
|
if (cGetMin < 1)
|
|
{
|
|
nextConsumer--;
|
|
continue;
|
|
}
|
|
int cGetMax = cellsPerNode - nextConsumer->has;
|
|
|
|
int p = nextProducer->procId;
|
|
|
|
int pSendMax = nextProducer->has - minCells;
|
|
|
|
if (pSendMax < 1)
|
|
{
|
|
nextProducer++;
|
|
continue;
|
|
}
|
|
|
|
int transferSize = (pSendMax < cGetMax) ? pSendMax : cGetMax;
|
|
|
|
if (me == p)
|
|
{
|
|
vtkIdType startCellId = nextProducer->had - nextProducer->has;
|
|
sendCells[c] = vtkIdList::New();
|
|
sendCells[c]->SetNumberOfIds(transferSize);
|
|
for (i=0; i<transferSize; i++)
|
|
{
|
|
sendCells[c]->SetId(i, startCellId++);
|
|
}
|
|
|
|
numTransferCells += transferSize;
|
|
}
|
|
|
|
nextProducer->has -= transferSize;
|
|
nextConsumer->has += transferSize;
|
|
|
|
continue;
|
|
}
|
|
|
|
delete [] procInfo;
|
|
|
|
if (sanityCheck > nprocsSquared)
|
|
{
|
|
vtkErrorMacro(<< "TestFixTooFewInputFiles error");
|
|
for (i=0; i<nprocs; i++)
|
|
{
|
|
if (sendCells[i])
|
|
{
|
|
sendCells[i]->Delete();
|
|
}
|
|
}
|
|
delete [] sendCells;
|
|
delete [] nodeType;
|
|
sendCells = NULL;
|
|
}
|
|
else if (nodeType[me] == Producer)
|
|
{
|
|
int keepCells = numMyCells - numTransferCells;
|
|
vtkIdType startCellId = (vtkIdType)numTransferCells;
|
|
sendCells[me] = vtkIdList::New();
|
|
sendCells[me]->SetNumberOfIds(keepCells);
|
|
for (i=0; i<keepCells; i++)
|
|
{
|
|
sendCells[me]->SetId(i, startCellId++);
|
|
}
|
|
}
|
|
}
|
|
|
|
vtkUnstructuredGrid *newGrid = NULL;
|
|
|
|
if (sendCells)
|
|
{
|
|
newGrid = this->ExchangeMergeSubGrids(
|
|
sendCells, DeleteYes, input, DeleteNo,
|
|
DuplicateCellsNo, GhostCellsNo, 0x0011);
|
|
|
|
delete [] sendCells;
|
|
delete [] nodeType;
|
|
}
|
|
|
|
return newGrid;
|
|
}
|
|
//-------------------------------------------------------------------------
|
|
// Communication routines - two versions:
|
|
// *Lean version use minimal memory
|
|
// *Fast versions use more memory, but are much faster
|
|
//-------------------------------------------------------------------------
|
|
|
|
void vtkDistributedDataFilter::SetUpPairWiseExchange()
|
|
{
|
|
int iam = this->MyId;
|
|
int nprocs = this->NumProcesses;
|
|
|
|
if (this->Target)
|
|
{
|
|
delete [] this->Target;
|
|
this->Target = NULL;
|
|
}
|
|
|
|
if (this->Source)
|
|
{
|
|
delete [] this->Source;
|
|
this->Source = NULL;
|
|
}
|
|
|
|
if (nprocs == 1)
|
|
{
|
|
return;
|
|
}
|
|
|
|
this->Target = new int [nprocs - 1];
|
|
this->Source = new int [nprocs - 1];
|
|
|
|
for (int i=1; i< nprocs; i++)
|
|
{
|
|
this->Target[i-1] = (iam + i) % nprocs;
|
|
this->Source[i-1] = (iam + nprocs - i) % nprocs;
|
|
}
|
|
}
|
|
void vtkDistributedDataFilter::FreeIntArrays(vtkIntArray **ar)
|
|
{
|
|
for (int i=0; i<this->NumProcesses; i++)
|
|
{
|
|
if (ar[i])
|
|
{
|
|
ar[i]->Delete();
|
|
}
|
|
}
|
|
|
|
delete [] ar;
|
|
}
|
|
void vtkDistributedDataFilter::FreeIdLists(vtkIdList**lists, int nlists)
|
|
{
|
|
for (int i=0; i<nlists; i++)
|
|
{
|
|
if (lists[i])
|
|
{
|
|
lists[i]->Delete();
|
|
lists[i] = NULL;
|
|
}
|
|
}
|
|
}
|
|
vtkIdType vtkDistributedDataFilter::GetIdListSize(vtkIdList **lists, int nlists)
|
|
{
|
|
vtkIdType numCells = 0;
|
|
|
|
for (int i=0; i<nlists; i++)
|
|
{
|
|
if (lists[i])
|
|
{
|
|
numCells += lists[i]->GetNumberOfIds();
|
|
}
|
|
}
|
|
|
|
return numCells;
|
|
}
|
|
vtkUnstructuredGrid *
|
|
vtkDistributedDataFilter::ExchangeMergeSubGrids(
|
|
vtkIdList **cellIds, int deleteCellIds,
|
|
vtkDataSet *myGrid, int deleteMyGrid,
|
|
int filterOutDuplicateCells, int ghostCellFlag,
|
|
int tag)
|
|
{
|
|
int nprocs = this->NumProcesses;
|
|
|
|
int *numLists = new int [nprocs];
|
|
|
|
vtkIdList ***listOfLists = new vtkIdList ** [nprocs];
|
|
|
|
for (int i=0; i<nprocs; i++)
|
|
{
|
|
if (cellIds[i] == NULL)
|
|
{
|
|
numLists[i] = 0;
|
|
}
|
|
else
|
|
{
|
|
numLists[i] = 1;
|
|
}
|
|
|
|
listOfLists[i] = &cellIds[i];
|
|
}
|
|
|
|
vtkUnstructuredGrid *grid = NULL;
|
|
|
|
if (this->UseMinimalMemory)
|
|
{
|
|
grid = this->ExchangeMergeSubGridsLean(listOfLists, numLists, deleteCellIds,
|
|
myGrid, deleteMyGrid, filterOutDuplicateCells, ghostCellFlag, tag);
|
|
}
|
|
else
|
|
{
|
|
grid = this->ExchangeMergeSubGridsFast(listOfLists, numLists, deleteCellIds,
|
|
myGrid, deleteMyGrid, filterOutDuplicateCells, ghostCellFlag, tag);
|
|
}
|
|
|
|
delete [] numLists;
|
|
delete [] listOfLists;
|
|
|
|
return grid;
|
|
}
|
|
vtkUnstructuredGrid *
|
|
vtkDistributedDataFilter::ExchangeMergeSubGrids(
|
|
vtkIdList ***cellIds, int *numLists, int deleteCellIds,
|
|
vtkDataSet *myGrid, int deleteMyGrid,
|
|
int filterOutDuplicateCells, int ghostCellFlag,
|
|
int tag)
|
|
{
|
|
vtkUnstructuredGrid *grid = NULL;
|
|
|
|
if (this->UseMinimalMemory)
|
|
{
|
|
grid = this->ExchangeMergeSubGridsLean(cellIds, numLists, deleteCellIds,
|
|
myGrid, deleteMyGrid, filterOutDuplicateCells, ghostCellFlag, tag);
|
|
}
|
|
else
|
|
{
|
|
grid = this->ExchangeMergeSubGridsFast(cellIds, numLists, deleteCellIds,
|
|
myGrid, deleteMyGrid, filterOutDuplicateCells, ghostCellFlag, tag);
|
|
}
|
|
return grid;
|
|
}
|
|
vtkIntArray *vtkDistributedDataFilter::ExchangeCounts(int myCount, int tag)
|
|
{
|
|
vtkIntArray *ia;
|
|
|
|
if (this->UseMinimalMemory)
|
|
{
|
|
ia = this->ExchangeCountsLean(myCount, tag);
|
|
}
|
|
else
|
|
{
|
|
ia = this->ExchangeCountsFast(myCount, tag);
|
|
}
|
|
return ia;
|
|
}
|
|
vtkFloatArray **vtkDistributedDataFilter::
|
|
ExchangeFloatArrays(vtkFloatArray **myArray, int deleteSendArrays, int tag)
|
|
{
|
|
vtkFloatArray **fa;
|
|
|
|
if (this->UseMinimalMemory)
|
|
{
|
|
fa = this->ExchangeFloatArraysLean(myArray, deleteSendArrays, tag);
|
|
}
|
|
else
|
|
{
|
|
fa = this->ExchangeFloatArraysFast(myArray, deleteSendArrays, tag);
|
|
}
|
|
return fa;
|
|
}
|
|
vtkIntArray **vtkDistributedDataFilter::
|
|
ExchangeIntArrays(vtkIntArray **myArray, int deleteSendArrays, int tag)
|
|
{
|
|
vtkIntArray **ia;
|
|
|
|
if (this->UseMinimalMemory)
|
|
{
|
|
ia = this->ExchangeIntArraysLean(myArray, deleteSendArrays, tag);
|
|
}
|
|
else
|
|
{
|
|
ia = this->ExchangeIntArraysFast(myArray, deleteSendArrays, tag);
|
|
}
|
|
return ia;
|
|
}
|
|
// ----------------------- Lean versions ----------------------------//
|
|
vtkIntArray *vtkDistributedDataFilter::ExchangeCountsLean(int myCount, int tag)
|
|
{
|
|
vtkIntArray *countArray = NULL;
|
|
|
|
#ifdef VTK_USE_MPI
|
|
int i;
|
|
int nprocs = this->NumProcesses;
|
|
|
|
vtkMPICommunicator::Request req;
|
|
vtkMPIController *mpiContr = vtkMPIController::SafeDownCast(this->Controller);
|
|
|
|
int *counts = new int [nprocs];
|
|
counts[this->MyId] = myCount;
|
|
|
|
if (!this->Source)
|
|
{
|
|
this->SetUpPairWiseExchange();
|
|
}
|
|
|
|
for (i = 0; i < this->NumProcesses - 1; i++)
|
|
{
|
|
int source = this->Source[i];
|
|
int target = this->Target[i];
|
|
mpiContr->NoBlockReceive(counts + source, 1, source, tag, req);
|
|
mpiContr->Send(&myCount, 1, target, tag);
|
|
req.Wait();
|
|
}
|
|
|
|
countArray = vtkIntArray::New();
|
|
countArray->SetArray(counts, nprocs, 0);
|
|
|
|
#else
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::ExchangeCounts requires MPI");
|
|
(void)myCount;
|
|
(void)tag;
|
|
#endif
|
|
|
|
return countArray;
|
|
}
|
|
vtkFloatArray **
|
|
vtkDistributedDataFilter::ExchangeFloatArraysLean(vtkFloatArray **myArray,
|
|
int deleteSendArrays, int tag)
|
|
{
|
|
vtkFloatArray **remoteArrays = NULL;
|
|
|
|
#ifdef VTK_USE_MPI
|
|
int i;
|
|
int nprocs = this->NumProcesses;
|
|
int me = this->MyId;
|
|
|
|
vtkMPICommunicator::Request req;
|
|
vtkMPIController *mpiContr = vtkMPIController::SafeDownCast(this->Controller);
|
|
|
|
int *recvSize = new int [nprocs];
|
|
int *sendSize = new int [nprocs];
|
|
|
|
if (!this->Source)
|
|
{
|
|
this->SetUpPairWiseExchange();
|
|
}
|
|
|
|
for (i= 0; i< nprocs; i++)
|
|
{
|
|
sendSize[i] = myArray[i] ? myArray[i]->GetNumberOfTuples() : 0;
|
|
recvSize[i] = 0;
|
|
}
|
|
|
|
// Exchange sizes
|
|
|
|
int nothers = nprocs - 1;
|
|
|
|
for (i = 0; i < nothers; i++)
|
|
{
|
|
int source = this->Source[i];
|
|
int target = this->Target[i];
|
|
mpiContr->NoBlockReceive(recvSize + source, 1, source, tag, req);
|
|
mpiContr->Send(sendSize + target, 1, target, tag);
|
|
req.Wait();
|
|
}
|
|
|
|
// Exchange int arrays
|
|
|
|
float **recvArrays = new float * [nprocs];
|
|
memset(recvArrays, 0, sizeof(float *) * nprocs);
|
|
|
|
if (sendSize[me] > 0) // sent myself an array
|
|
{
|
|
recvSize[me] = sendSize[me];
|
|
recvArrays[me] = new float [sendSize[me]];
|
|
memcpy(recvArrays[me], myArray[me]->GetPointer(0), sendSize[me] * sizeof(float));
|
|
}
|
|
|
|
for (i = 0; i < nothers; i++)
|
|
{
|
|
int source = this->Source[i];
|
|
int target = this->Target[i];
|
|
recvArrays[source] = NULL;
|
|
|
|
if (recvSize[source] > 0)
|
|
{
|
|
recvArrays[source] = new float [recvSize[source]];
|
|
if (recvArrays[source] == NULL)
|
|
{
|
|
vtkErrorMacro(<<
|
|
"vtkDistributedDataFilter::ExchangeIntArrays memory allocation");
|
|
return NULL;
|
|
}
|
|
mpiContr->NoBlockReceive(recvArrays[source], recvSize[source], source, tag, req);
|
|
}
|
|
|
|
if (sendSize[target] > 0)
|
|
{
|
|
mpiContr->Send(myArray[target]->GetPointer(0), sendSize[target], target, tag);
|
|
}
|
|
|
|
if (myArray[target] && deleteSendArrays)
|
|
{
|
|
myArray[target]->Delete();
|
|
}
|
|
|
|
if (recvSize[source] > 0)
|
|
{
|
|
req.Wait();
|
|
}
|
|
}
|
|
|
|
if (deleteSendArrays)
|
|
{
|
|
if (myArray[me])
|
|
{
|
|
myArray[me]->Delete();
|
|
}
|
|
delete [] myArray;
|
|
}
|
|
|
|
delete [] sendSize;
|
|
|
|
remoteArrays = new vtkFloatArray * [nprocs];
|
|
|
|
for (i=0; i<nprocs; i++)
|
|
{
|
|
if (recvSize[i] > 0)
|
|
{
|
|
remoteArrays[i] = vtkFloatArray::New();
|
|
remoteArrays[i]->SetArray(recvArrays[i], recvSize[i], 0);
|
|
}
|
|
else
|
|
{
|
|
remoteArrays[i] = NULL;
|
|
}
|
|
}
|
|
|
|
delete [] recvArrays;
|
|
delete [] recvSize;
|
|
|
|
#else
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::ExchangeFloatArrays requires MPI");
|
|
(void)myArray;
|
|
(void)deleteSendArrays;
|
|
(void)tag;
|
|
#endif
|
|
|
|
return remoteArrays;
|
|
}
|
|
vtkIntArray **
|
|
vtkDistributedDataFilter::ExchangeIntArraysLean(vtkIntArray **myArray,
|
|
int deleteSendArrays, int tag)
|
|
{
|
|
vtkIntArray **remoteArrays = NULL;
|
|
|
|
#ifdef VTK_USE_MPI
|
|
int i;
|
|
int nprocs = this->NumProcesses;
|
|
int me = this->MyId;
|
|
|
|
vtkMPICommunicator::Request req;
|
|
vtkMPIController *mpiContr = vtkMPIController::SafeDownCast(this->Controller);
|
|
|
|
int *recvSize = new int [nprocs];
|
|
int *sendSize = new int [nprocs];
|
|
|
|
if (!this->Source)
|
|
{
|
|
this->SetUpPairWiseExchange();
|
|
}
|
|
|
|
for (i= 0; i< nprocs; i++)
|
|
{
|
|
sendSize[i] = myArray[i] ? myArray[i]->GetNumberOfTuples() : 0;
|
|
recvSize[i] = 0;
|
|
}
|
|
|
|
// Exchange sizes
|
|
|
|
int nothers = nprocs - 1;
|
|
|
|
for (i = 0; i < nothers; i++)
|
|
{
|
|
int source = this->Source[i];
|
|
int target = this->Target[i];
|
|
mpiContr->NoBlockReceive(recvSize + source, 1, source, tag, req);
|
|
mpiContr->Send(sendSize + target, 1, target, tag);
|
|
req.Wait();
|
|
}
|
|
|
|
// Exchange int arrays
|
|
|
|
int **recvArrays = new int * [nprocs];
|
|
memset(recvArrays, 0, sizeof(int *) * nprocs);
|
|
|
|
if (sendSize[me] > 0) // sent myself an array
|
|
{
|
|
recvSize[me] = sendSize[me];
|
|
recvArrays[me] = new int [sendSize[me]];
|
|
memcpy(recvArrays[me], myArray[me]->GetPointer(0), sendSize[me] * sizeof(int));
|
|
}
|
|
|
|
for (i = 0; i < nothers; i++)
|
|
{
|
|
int source = this->Source[i];
|
|
int target = this->Target[i];
|
|
recvArrays[source] = NULL;
|
|
|
|
if (recvSize[source] > 0)
|
|
{
|
|
recvArrays[source] = new int [recvSize[source]];
|
|
if (recvArrays[source] == NULL)
|
|
{
|
|
vtkErrorMacro(<<
|
|
"vtkDistributedDataFilter::ExchangeIntArrays memory allocation");
|
|
return NULL;
|
|
}
|
|
mpiContr->NoBlockReceive(recvArrays[source], recvSize[source], source, tag, req);
|
|
}
|
|
|
|
if (sendSize[target] > 0)
|
|
{
|
|
mpiContr->Send(myArray[target]->GetPointer(0), sendSize[target], target, tag);
|
|
}
|
|
|
|
if (myArray[target] && deleteSendArrays)
|
|
{
|
|
myArray[target]->Delete();
|
|
}
|
|
|
|
if (recvSize[source] > 0)
|
|
{
|
|
req.Wait();
|
|
}
|
|
}
|
|
|
|
if (deleteSendArrays)
|
|
{
|
|
if (myArray[me])
|
|
{
|
|
myArray[me]->Delete();
|
|
}
|
|
delete [] myArray;
|
|
}
|
|
|
|
delete [] sendSize;
|
|
|
|
remoteArrays = new vtkIntArray * [nprocs];
|
|
|
|
for (i=0; i<nprocs; i++)
|
|
{
|
|
if (recvSize[i] > 0)
|
|
{
|
|
remoteArrays[i] = vtkIntArray::New();
|
|
remoteArrays[i]->SetArray(recvArrays[i], recvSize[i], 0);
|
|
}
|
|
else
|
|
{
|
|
remoteArrays[i] = NULL;
|
|
}
|
|
}
|
|
|
|
delete [] recvArrays;
|
|
delete [] recvSize;
|
|
|
|
#else
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::ExchangeIntArrays requires MPI");
|
|
(void)myArray;
|
|
(void)deleteSendArrays;
|
|
(void)tag;
|
|
#endif
|
|
|
|
return remoteArrays;
|
|
}
|
|
vtkUnstructuredGrid *
|
|
vtkDistributedDataFilter::ExchangeMergeSubGridsLean(
|
|
vtkIdList ***cellIds, int *numLists, int deleteCellIds,
|
|
vtkDataSet *myGrid, int deleteMyGrid,
|
|
int filterOutDuplicateCells, // flag if different processes may send same cells
|
|
int ghostCellFlag, // flag if these cells are ghost cells
|
|
int tag)
|
|
{
|
|
vtkUnstructuredGrid *mergedGrid = NULL;
|
|
#ifdef VTK_USE_MPI
|
|
int i;
|
|
int packedGridSendSize=0, packedGridRecvSize=0;
|
|
char *packedGridSend=NULL, *packedGridRecv=NULL;
|
|
int recvBufSize=0;
|
|
int numReceivedGrids = 0;
|
|
|
|
int nprocs = this->NumProcesses;
|
|
int iam = this->MyId;
|
|
|
|
vtkMPIController *mpiContr = vtkMPIController::SafeDownCast(this->Controller);
|
|
vtkMPICommunicator::Request req;
|
|
|
|
vtkDataSet *tmpGrid = myGrid->NewInstance();
|
|
tmpGrid->ShallowCopy(myGrid);
|
|
|
|
vtkModelMetadata *mmd = NULL;
|
|
|
|
if (vtkDistributedDataFilter::HasMetadata(myGrid) && !ghostCellFlag)
|
|
{
|
|
// Pull metadata out of grid
|
|
|
|
mmd = vtkModelMetadata::New();
|
|
mmd->Unpack(tmpGrid, DeleteYes);
|
|
}
|
|
|
|
vtkDataSet **grids = new vtkDataSet * [nprocs];
|
|
|
|
if (numLists[iam] > 0)
|
|
{
|
|
// I was extracting/packing/sending/unpacking ugrids of zero cells,
|
|
// and this caused corrupted data structures. I don't know why, but
|
|
// I am now being careful not to do that.
|
|
|
|
vtkIdType numCells =
|
|
vtkDistributedDataFilter::GetIdListSize(cellIds[iam], numLists[iam]);
|
|
|
|
if (numCells > 0)
|
|
{
|
|
grids[numReceivedGrids++] =
|
|
this->ExtractCells(cellIds[iam], numLists[iam], deleteCellIds, tmpGrid, mmd);
|
|
}
|
|
else if (deleteCellIds)
|
|
{
|
|
vtkDistributedDataFilter::FreeIdLists(cellIds[iam], numLists[iam]);
|
|
}
|
|
}
|
|
|
|
if (this->Source == NULL)
|
|
{
|
|
this->SetUpPairWiseExchange();
|
|
}
|
|
|
|
int nothers = nprocs - 1;
|
|
|
|
for (i=0; i<nothers; i++)
|
|
{
|
|
int target = this->Target[i];
|
|
int source = this->Source[i];
|
|
|
|
packedGridSendSize = 0;
|
|
|
|
if (cellIds[target] && (numLists[target] > 0))
|
|
{
|
|
vtkIdType numCells = vtkDistributedDataFilter::GetIdListSize(
|
|
cellIds[target], numLists[target]);
|
|
|
|
if (numCells > 0)
|
|
{
|
|
vtkUnstructuredGrid *sendGrid =
|
|
this->ExtractCells(cellIds[target], numLists[target],
|
|
deleteCellIds, tmpGrid, mmd);
|
|
|
|
packedGridSend = this->MarshallDataSet(sendGrid, packedGridSendSize);
|
|
sendGrid->Delete();
|
|
}
|
|
else if (deleteCellIds)
|
|
{
|
|
vtkDistributedDataFilter::FreeIdLists(cellIds[target], numLists[target]);
|
|
}
|
|
}
|
|
|
|
// exchange size of packed grids
|
|
|
|
mpiContr->NoBlockReceive(&packedGridRecvSize, 1, source, tag, req);
|
|
mpiContr->Send(&packedGridSendSize, 1, target, tag);
|
|
req.Wait();
|
|
|
|
if (packedGridRecvSize > recvBufSize)
|
|
{
|
|
if (packedGridRecv)
|
|
{
|
|
delete [] packedGridRecv;
|
|
}
|
|
packedGridRecv = new char [packedGridRecvSize];
|
|
if (!packedGridRecv)
|
|
{
|
|
vtkErrorMacro(<<
|
|
"vtkDistributedDataFilter::ExchangeMergeSubGrids memory allocation");
|
|
return NULL;
|
|
}
|
|
recvBufSize = packedGridRecvSize;
|
|
}
|
|
|
|
if (packedGridRecvSize > 0)
|
|
{
|
|
mpiContr->NoBlockReceive(packedGridRecv, packedGridRecvSize, source, tag, req);
|
|
}
|
|
|
|
if (packedGridSendSize > 0)
|
|
{
|
|
mpiContr->Send(packedGridSend, packedGridSendSize, target, tag);
|
|
delete [] packedGridSend;
|
|
}
|
|
|
|
if (packedGridRecvSize > 0)
|
|
{
|
|
req.Wait();
|
|
|
|
grids[numReceivedGrids++] =
|
|
this->UnMarshallDataSet(packedGridRecv, packedGridRecvSize);
|
|
}
|
|
}
|
|
|
|
tmpGrid->Delete();
|
|
|
|
if (recvBufSize > 0)
|
|
{
|
|
delete [] packedGridRecv;
|
|
packedGridRecv = NULL;
|
|
}
|
|
|
|
if (numReceivedGrids > 1)
|
|
{
|
|
// Merge received grids
|
|
|
|
const char *globalNodeIds = this->GetGlobalNodeIdArrayName(myGrid);
|
|
const char *globalElementIds = NULL;
|
|
|
|
if (filterOutDuplicateCells)
|
|
{
|
|
globalElementIds = this->GetGlobalElementIdArrayName(myGrid);
|
|
}
|
|
|
|
// this call will merge the grids and then delete them
|
|
|
|
float tolerance = 0.0;
|
|
|
|
if (this->Kdtree)
|
|
{
|
|
tolerance = (float)this->Kdtree->GetFudgeFactor();
|
|
}
|
|
|
|
mergedGrid =
|
|
vtkDistributedDataFilter::MergeGrids(grids, numReceivedGrids, DeleteYes,
|
|
globalNodeIds, tolerance, globalElementIds);
|
|
|
|
}
|
|
else if (numReceivedGrids == 1)
|
|
{
|
|
mergedGrid = vtkUnstructuredGrid::SafeDownCast(grids[0]);
|
|
}
|
|
else
|
|
{
|
|
mergedGrid = this->ExtractZeroCellGrid(myGrid, mmd);
|
|
}
|
|
|
|
if (mmd)
|
|
{
|
|
mmd->Delete();
|
|
}
|
|
|
|
if (deleteMyGrid)
|
|
{
|
|
myGrid->Delete();
|
|
}
|
|
|
|
delete [] grids;
|
|
|
|
#else
|
|
(void)cellIds; // This is just here for successful compilation,
|
|
(void)numLists; // it will never execute. If !VTK_USE_MPI, we
|
|
(void)deleteCellIds; // never get this far.
|
|
(void)myGrid;
|
|
(void)deleteMyGrid;
|
|
(void)filterOutDuplicateCells;
|
|
(void)tag;
|
|
(void)ghostCellFlag;
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::ExchangeMergeSubGrids requires MPI");
|
|
#endif
|
|
|
|
return mergedGrid;
|
|
}
|
|
// ----------------------- Fast versions ----------------------------//
|
|
vtkIntArray *vtkDistributedDataFilter::ExchangeCountsFast(int myCount, int tag)
|
|
{
|
|
vtkIntArray *countArray = NULL;
|
|
|
|
#ifdef VTK_USE_MPI
|
|
int i;
|
|
int nprocs = this->NumProcesses;
|
|
int me = this->MyId;
|
|
|
|
vtkMPICommunicator::Request *req = new vtkMPICommunicator::Request [nprocs];
|
|
vtkMPIController *mpiContr = vtkMPIController::SafeDownCast(this->Controller);
|
|
|
|
int *counts = new int [nprocs];
|
|
counts[me] = myCount;
|
|
|
|
for (i = 0; i < nprocs; i++)
|
|
{
|
|
if (i == me)
|
|
{
|
|
continue;
|
|
}
|
|
mpiContr->NoBlockReceive(counts + i, 1, i, tag, req[i]);
|
|
}
|
|
|
|
mpiContr->Barrier();
|
|
|
|
for (i = 0; i < nprocs; i++)
|
|
{
|
|
if (i == me)
|
|
{
|
|
continue;
|
|
}
|
|
mpiContr->Send(&myCount, 1, i, tag);
|
|
}
|
|
|
|
countArray = vtkIntArray::New();
|
|
countArray->SetArray(counts, nprocs, 0);
|
|
|
|
for (i = 0; i < nprocs; i++)
|
|
{
|
|
if (i == me)
|
|
{
|
|
continue;
|
|
}
|
|
req[i].Wait();
|
|
}
|
|
|
|
delete [] req;
|
|
|
|
#else
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::ExchangeCounts requires MPI");
|
|
(void)myCount;
|
|
(void)tag;
|
|
#endif
|
|
|
|
return countArray;
|
|
}
|
|
vtkFloatArray **
|
|
vtkDistributedDataFilter::ExchangeFloatArraysFast(vtkFloatArray **myArray,
|
|
int deleteSendArrays, int tag)
|
|
{
|
|
vtkFloatArray **fa = NULL;
|
|
#ifdef VTK_USE_MPI
|
|
int proc;
|
|
int nprocs = this->NumProcesses;
|
|
int iam = this->MyId;
|
|
|
|
vtkMPIController *mpiContr = vtkMPIController::SafeDownCast(this->Controller);
|
|
|
|
int *sendSize = new int [nprocs];
|
|
int *recvSize = new int [nprocs];
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
recvSize[proc] = sendSize[proc] = 0;
|
|
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (myArray[proc])
|
|
{
|
|
sendSize[proc] = myArray[proc]->GetNumberOfTuples();
|
|
}
|
|
}
|
|
|
|
// Exchange sizes of arrays to send and receive
|
|
|
|
vtkMPICommunicator::Request *reqBuf = new vtkMPICommunicator::Request [nprocs];
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
mpiContr->NoBlockReceive(recvSize + proc, 1, proc, tag, reqBuf[proc]);
|
|
}
|
|
|
|
mpiContr->Barrier();
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
mpiContr->Send(sendSize + proc, 1, proc, tag);
|
|
}
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
reqBuf[proc].Wait();
|
|
}
|
|
|
|
// Allocate buffers and post receives
|
|
|
|
float **recvBufs = new float * [nprocs];
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (recvSize[proc] > 0)
|
|
{
|
|
recvBufs[proc] = new float [recvSize[proc]];
|
|
mpiContr->NoBlockReceive(recvBufs[proc], recvSize[proc], proc, tag, reqBuf[proc]);
|
|
}
|
|
else
|
|
{
|
|
recvBufs[proc] = NULL;
|
|
}
|
|
}
|
|
|
|
mpiContr->Barrier();
|
|
|
|
// Send all arrays
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (sendSize[proc] > 0)
|
|
{
|
|
mpiContr->Send(myArray[proc]->GetPointer(0), sendSize[proc], proc, tag);
|
|
}
|
|
}
|
|
delete [] sendSize;
|
|
|
|
// If I want to send an array to myself, place it in output now
|
|
|
|
if (myArray[iam])
|
|
{
|
|
recvSize[iam] = myArray[iam]->GetNumberOfTuples();
|
|
if (recvSize[iam] > 0)
|
|
{
|
|
recvBufs[iam] = new float [recvSize[iam]];
|
|
memcpy(recvBufs[iam], myArray[iam]->GetPointer(0), recvSize[iam] * sizeof(float));
|
|
}
|
|
}
|
|
|
|
if (deleteSendArrays)
|
|
{
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (myArray[proc])
|
|
{
|
|
myArray[proc]->Delete();
|
|
}
|
|
}
|
|
delete [] myArray;
|
|
}
|
|
|
|
// Await incoming arrays
|
|
|
|
fa = new vtkFloatArray * [nprocs];
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (recvBufs[proc])
|
|
{
|
|
fa[proc] = vtkFloatArray::New();
|
|
fa[proc]->SetArray(recvBufs[proc], recvSize[proc], 0);
|
|
}
|
|
else
|
|
{
|
|
fa[proc] = NULL;
|
|
}
|
|
}
|
|
|
|
delete [] recvSize;
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
if (recvBufs[proc])
|
|
{
|
|
reqBuf[proc].Wait();
|
|
}
|
|
}
|
|
|
|
delete [] reqBuf;
|
|
delete [] recvBufs;
|
|
|
|
#else
|
|
(void)myArray;
|
|
(void)deleteSendArrays;
|
|
(void)tag;
|
|
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::ExchangeFloatArrays requires MPI");
|
|
#endif
|
|
|
|
return fa;
|
|
}
|
|
vtkIntArray **
|
|
vtkDistributedDataFilter::ExchangeIntArraysFast(vtkIntArray **myArray,
|
|
int deleteSendArrays, int tag)
|
|
{
|
|
vtkIntArray **ia = NULL;
|
|
#ifdef VTK_USE_MPI
|
|
int proc;
|
|
int nprocs = this->NumProcesses;
|
|
int iam = this->MyId;
|
|
|
|
vtkMPIController *mpiContr = vtkMPIController::SafeDownCast(this->Controller);
|
|
|
|
int *sendSize = new int [nprocs];
|
|
int *recvSize = new int [nprocs];
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
recvSize[proc] = sendSize[proc] = 0;
|
|
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (myArray[proc])
|
|
{
|
|
sendSize[proc] = myArray[proc]->GetNumberOfTuples();
|
|
}
|
|
}
|
|
|
|
// Exchange sizes of arrays to send and receive
|
|
|
|
vtkMPICommunicator::Request *reqBuf = new vtkMPICommunicator::Request [nprocs];
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
mpiContr->NoBlockReceive(recvSize + proc, 1, proc, tag, reqBuf[proc]);
|
|
}
|
|
|
|
mpiContr->Barrier();
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
mpiContr->Send(sendSize + proc, 1, proc, tag);
|
|
}
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
reqBuf[proc].Wait();
|
|
}
|
|
|
|
// Allocate buffers and post receives
|
|
|
|
int **recvBufs = new int * [nprocs];
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (recvSize[proc] > 0)
|
|
{
|
|
recvBufs[proc] = new int [recvSize[proc]];
|
|
mpiContr->NoBlockReceive(recvBufs[proc], recvSize[proc], proc, tag, reqBuf[proc]);
|
|
}
|
|
else
|
|
{
|
|
recvBufs[proc] = NULL;
|
|
}
|
|
}
|
|
|
|
mpiContr->Barrier();
|
|
|
|
// Send all arrays
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (sendSize[proc] > 0)
|
|
{
|
|
mpiContr->Send(myArray[proc]->GetPointer(0), sendSize[proc], proc, tag);
|
|
}
|
|
}
|
|
delete [] sendSize;
|
|
|
|
// If I want to send an array to myself, place it in output now
|
|
|
|
if (myArray[iam])
|
|
{
|
|
recvSize[iam] = myArray[iam]->GetNumberOfTuples();
|
|
if (recvSize[iam] > 0)
|
|
{
|
|
recvBufs[iam] = new int [recvSize[iam]];
|
|
memcpy(recvBufs[iam], myArray[iam]->GetPointer(0), recvSize[iam] * sizeof(int));
|
|
}
|
|
}
|
|
|
|
if (deleteSendArrays)
|
|
{
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (myArray[proc])
|
|
{
|
|
myArray[proc]->Delete();
|
|
}
|
|
}
|
|
delete [] myArray;
|
|
}
|
|
|
|
// Await incoming arrays
|
|
|
|
ia = new vtkIntArray * [nprocs];
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (recvBufs[proc])
|
|
{
|
|
ia[proc] = vtkIntArray::New();
|
|
ia[proc]->SetArray(recvBufs[proc], recvSize[proc], 0);
|
|
}
|
|
else
|
|
{
|
|
ia[proc] = NULL;
|
|
}
|
|
}
|
|
|
|
delete [] recvSize;
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
if (recvBufs[proc])
|
|
{
|
|
reqBuf[proc].Wait();
|
|
}
|
|
}
|
|
|
|
delete [] reqBuf;
|
|
delete [] recvBufs;
|
|
|
|
#else
|
|
(void)myArray;
|
|
(void)deleteSendArrays;
|
|
(void)tag;
|
|
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::ExchangeIntArrays requires MPI");
|
|
#endif
|
|
|
|
return ia;
|
|
}
|
|
vtkUnstructuredGrid *
|
|
vtkDistributedDataFilter::ExchangeMergeSubGridsFast(
|
|
vtkIdList ***cellIds, int *numLists, int deleteCellIds,
|
|
vtkDataSet *myGrid, int deleteMyGrid,
|
|
int filterOutDuplicateCells, // flag if different processes may send same cells
|
|
int ghostCellFlag, // flag if these are ghost cells
|
|
int tag)
|
|
{
|
|
vtkUnstructuredGrid *mergedGrid = NULL;
|
|
#ifdef VTK_USE_MPI
|
|
int proc;
|
|
int nprocs = this->NumProcesses;
|
|
int iam = this->MyId;
|
|
|
|
vtkMPIController *mpiContr = vtkMPIController::SafeDownCast(this->Controller);
|
|
|
|
vtkUnstructuredGrid **grids = new vtkUnstructuredGrid * [nprocs];
|
|
char **sendBufs = new char * [nprocs];
|
|
char **recvBufs = new char * [nprocs];
|
|
int *sendSize = new int [nprocs];
|
|
int *recvSize = new int [nprocs];
|
|
|
|
// create & pack all sub grids
|
|
|
|
vtkDataSet *tmpGrid = myGrid->NewInstance();
|
|
tmpGrid->ShallowCopy(myGrid);
|
|
|
|
vtkModelMetadata *mmd = NULL;
|
|
|
|
if (vtkDistributedDataFilter::HasMetadata(tmpGrid) && !ghostCellFlag)
|
|
{
|
|
// Pull metadata out of grid
|
|
|
|
mmd = vtkModelMetadata::New();
|
|
mmd->Unpack(tmpGrid, DeleteYes);
|
|
}
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
recvSize[proc] = sendSize[proc] = 0;
|
|
grids[proc] = NULL;
|
|
sendBufs[proc] = recvBufs[proc] = NULL;
|
|
|
|
if (numLists[proc] > 0)
|
|
{
|
|
vtkIdType numCells =
|
|
vtkDistributedDataFilter::GetIdListSize(cellIds[proc], numLists[proc]);
|
|
|
|
if (numCells > 0)
|
|
{
|
|
grids[proc] =
|
|
vtkDistributedDataFilter::ExtractCells(cellIds[proc], numLists[proc],
|
|
deleteCellIds, tmpGrid, mmd);
|
|
|
|
if (proc != iam)
|
|
{
|
|
sendBufs[proc] = this->MarshallDataSet(grids[proc], sendSize[proc]);
|
|
grids[proc]->Delete();
|
|
grids[proc] = NULL;
|
|
}
|
|
}
|
|
else if (deleteCellIds)
|
|
{
|
|
vtkDistributedDataFilter::FreeIdLists(cellIds[proc], numLists[proc]);
|
|
}
|
|
}
|
|
}
|
|
|
|
tmpGrid->Delete();
|
|
|
|
// Exchange sizes of grids to send and receive
|
|
|
|
vtkMPICommunicator::Request *reqBuf = new vtkMPICommunicator::Request [nprocs];
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
mpiContr->NoBlockReceive(recvSize + proc, 1, proc, tag, reqBuf[proc]);
|
|
}
|
|
|
|
mpiContr->Barrier();
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
mpiContr->Send(sendSize + proc, 1, proc, tag);
|
|
}
|
|
|
|
for (proc=0; proc<nprocs; proc++)
|
|
{
|
|
if (proc == iam)
|
|
{
|
|
continue;
|
|
}
|
|
reqBuf[proc].Wait();
|
|
}
|
|
|
|
// Allocate buffers and post receives
|
|
|
|
int numReceives = 0;
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (recvSize[proc] > 0)
|
|
{
|
|
recvBufs[proc] = new char [recvSize[proc]];
|
|
mpiContr->NoBlockReceive(recvBufs[proc], recvSize[proc], proc, tag, reqBuf[proc]);
|
|
numReceives++;
|
|
}
|
|
}
|
|
|
|
mpiContr->Barrier();
|
|
|
|
// Send all sub grids, then delete them
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (sendSize[proc] > 0)
|
|
{
|
|
mpiContr->Send(sendBufs[proc], sendSize[proc], proc, tag);
|
|
}
|
|
}
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (sendSize[proc] > 0)
|
|
{
|
|
delete [] sendBufs[proc];
|
|
}
|
|
}
|
|
|
|
delete [] sendSize;
|
|
delete [] sendBufs;
|
|
|
|
// Await incoming sub grids, unpack them
|
|
|
|
while (numReceives > 0)
|
|
{
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (recvBufs[proc] && (reqBuf[proc].Test() == 1))
|
|
{
|
|
grids[proc] = this->UnMarshallDataSet(recvBufs[proc], recvSize[proc]);
|
|
delete [] recvBufs[proc];
|
|
recvBufs[proc] = NULL;
|
|
numReceives--;
|
|
}
|
|
}
|
|
}
|
|
|
|
delete [] reqBuf;
|
|
delete [] recvBufs;
|
|
delete [] recvSize;
|
|
|
|
// Merge received grids
|
|
|
|
float tolerance = 0.0;
|
|
|
|
if (this->Kdtree)
|
|
{
|
|
tolerance = (float)this->Kdtree->GetFudgeFactor();
|
|
}
|
|
|
|
int numReceivedGrids = 0;
|
|
|
|
vtkDataSet **ds = new vtkDataSet * [nprocs];
|
|
|
|
for (proc=0; proc < nprocs; proc++)
|
|
{
|
|
if (grids[proc] != NULL)
|
|
{
|
|
ds[numReceivedGrids++] = static_cast<vtkDataSet *>(grids[proc]);
|
|
}
|
|
}
|
|
|
|
delete [] grids;
|
|
|
|
if (numReceivedGrids > 1)
|
|
{
|
|
const char *globalNodeIds = this->GetGlobalNodeIdArrayName(ds[0]);
|
|
const char *globalCellIds = NULL;
|
|
if (filterOutDuplicateCells)
|
|
{
|
|
globalCellIds = this->GetGlobalElementIdArrayName(ds[0]);
|
|
}
|
|
|
|
// this call will merge the grids and then delete them
|
|
mergedGrid =
|
|
vtkDistributedDataFilter::MergeGrids(ds, numReceivedGrids, DeleteYes,
|
|
globalNodeIds, tolerance, globalCellIds);
|
|
|
|
}
|
|
else if (numReceivedGrids == 1)
|
|
{
|
|
mergedGrid = vtkUnstructuredGrid::SafeDownCast(ds[0]);
|
|
}
|
|
else
|
|
{
|
|
mergedGrid = this->ExtractZeroCellGrid(myGrid, mmd);
|
|
}
|
|
|
|
if (mmd)
|
|
{
|
|
mmd->Delete();
|
|
}
|
|
|
|
if (deleteMyGrid)
|
|
{
|
|
myGrid->Delete();
|
|
}
|
|
|
|
delete [] ds;
|
|
|
|
#else
|
|
(void)cellIds; // This is just here for successful compilation,
|
|
(void)numLists; // it will never execute. If !VTK_USE_MPI, we
|
|
(void)deleteCellIds; // never get this far.
|
|
(void)myGrid;
|
|
(void)deleteMyGrid;
|
|
(void)filterOutDuplicateCells;
|
|
(void)tag;
|
|
(void)ghostCellFlag;
|
|
|
|
vtkErrorMacro(<< "vtkDistributedDataFilter::ExchangeMergeSubGrids requires MPI");
|
|
#endif
|
|
|
|
return mergedGrid;
|
|
}
|
|
void vtkDistributedDataFilter::AddMetadata(vtkUnstructuredGrid *grid,
|
|
vtkModelMetadata *mmd)
|
|
{
|
|
const char *eltIdName = this->GetGlobalElementIdArrayName(grid);
|
|
vtkDataArray *da = grid->GetCellData()->GetArray(eltIdName);
|
|
vtkIntArray *ia = vtkIntArray::SafeDownCast(da);
|
|
|
|
const char *nodeIdName = this->GetGlobalNodeIdArrayName(grid);
|
|
|
|
// Extract the metadata for all cells in this grid
|
|
|
|
vtkModelMetadata *submmd =
|
|
mmd->ExtractModelMetadata(ia, // extract metadata for these cells
|
|
grid, // in this grid
|
|
eltIdName, // global cell ID array name
|
|
nodeIdName); // global node ID array name
|
|
|
|
// Pack that metadata into field arrays of the grid
|
|
|
|
submmd->Pack(grid);
|
|
submmd->Delete();
|
|
}
|
|
|
|
vtkUnstructuredGrid *vtkDistributedDataFilter::MPIRedistribute(vtkDataSet *in,
|
|
vtkDataSet *input)
|
|
{
|
|
int proc;
|
|
int nprocs = this->NumProcesses;
|
|
|
|
// A cell belongs to a spatial region if it's centroid lies in that
|
|
// region. The kdtree object can create a list for each region of the
|
|
// IDs of each cell I have read in that belong in that region. If we
|
|
// are building subgrids of all cells that intersect a region (a
|
|
// superset of all cells that belong to a region) then the kdtree object
|
|
// can build another set of lists of all cells that intersect each
|
|
// region (but don't have their centroid in that region).
|
|
|
|
if (this->IncludeAllIntersectingCells)
|
|
{
|
|
// TO DO:
|
|
// We actually compute whether a cell intersects a spatial region.
|
|
// This can be a lengthy calculation. Perhaps it's good enough
|
|
// to compute whether a cell's bounding box intersects the region.
|
|
// Some of the cells we list will actually not be in the region, but
|
|
// if we are clipping later, it doesn't matter.
|
|
//
|
|
// Is there any rendering algorithm that needs exactly all cells
|
|
// which intersect the region, and no more?
|
|
|
|
this->Kdtree->IncludeRegionBoundaryCellsOn(); // SLOW!!
|
|
}
|
|
|
|
this->Kdtree->CreateCellLists(); // required by GetCellIdsForProcess
|
|
|
|
vtkIdList ***procCellLists = new vtkIdList ** [nprocs];
|
|
int *numLists = new int [nprocs];
|
|
|
|
for (proc = 0; proc < this->NumProcesses; proc++)
|
|
{
|
|
procCellLists[proc] = this->GetCellIdsForProcess(proc, numLists + proc);
|
|
}
|
|
|
|
int deleteDataSet = DeleteNo;
|
|
|
|
if (in != input)
|
|
{
|
|
deleteDataSet = DeleteYes;
|
|
}
|
|
|
|
vtkUnstructuredGrid *myNewGrid =
|
|
this->ExchangeMergeSubGrids(procCellLists, numLists, DeleteNo,
|
|
in, deleteDataSet, DuplicateCellsNo, GhostCellsNo, 0x0012);
|
|
|
|
for (proc = 0; proc < nprocs; proc++)
|
|
{
|
|
delete [] procCellLists[proc];
|
|
}
|
|
|
|
delete [] procCellLists;
|
|
delete [] numLists;
|
|
|
|
if (myNewGrid && (this->GhostLevel > 0))
|
|
{
|
|
vtkDistributedDataFilter::AddConstantUnsignedCharCellArray(
|
|
myNewGrid, "vtkGhostLevels", 0);
|
|
vtkDistributedDataFilter::AddConstantUnsignedCharPointArray(
|
|
myNewGrid, "vtkGhostLevels", 0);
|
|
}
|
|
return myNewGrid;
|
|
}
|
|
|
|
char *vtkDistributedDataFilter::MarshallDataSet(vtkUnstructuredGrid *extractedGrid, int &len)
|
|
{
|
|
// taken from vtkCommunicator::WriteDataSet
|
|
|
|
vtkUnstructuredGrid *copy;
|
|
vtkDataSetWriter *writer = vtkDataSetWriter::New();
|
|
|
|
copy = extractedGrid->NewInstance();
|
|
copy->ShallowCopy(extractedGrid);
|
|
|
|
// There is a problem with binary files with no data.
|
|
if (copy->GetNumberOfCells() > 0)
|
|
{
|
|
writer->SetFileTypeToBinary();
|
|
}
|
|
writer->WriteToOutputStringOn();
|
|
writer->SetInput(copy);
|
|
|
|
writer->Write();
|
|
|
|
len = writer->GetOutputStringLength();
|
|
|
|
char *packedFormat = writer->RegisterAndGetOutputString();
|
|
|
|
writer->Delete();
|
|
|
|
copy->Delete();
|
|
|
|
return packedFormat;
|
|
}
|
|
vtkUnstructuredGrid *vtkDistributedDataFilter::UnMarshallDataSet(char *buf, int size)
|
|
{
|
|
// taken from vtkCommunicator::ReadDataSet
|
|
|
|
vtkDataSetReader *reader = vtkDataSetReader::New();
|
|
|
|
reader->ReadFromInputStringOn();
|
|
|
|
vtkCharArray* mystring = vtkCharArray::New();
|
|
|
|
mystring->SetArray(buf, size, 1);
|
|
|
|
reader->SetInputArray(mystring);
|
|
mystring->Delete();
|
|
|
|
vtkDataSet *output = reader->GetOutput();
|
|
output->Update();
|
|
|
|
vtkUnstructuredGrid *newGrid = vtkUnstructuredGrid::New();
|
|
|
|
newGrid->ShallowCopy(output);
|
|
|
|
reader->Delete();
|
|
|
|
return newGrid;
|
|
}
|
|
vtkUnstructuredGrid
|
|
*vtkDistributedDataFilter::ExtractCells(vtkIdList *cells, int deleteCellLists,
|
|
vtkDataSet *in, vtkModelMetadata *mmd)
|
|
{
|
|
vtkIdList *tempCellList = NULL;
|
|
|
|
if (cells == NULL)
|
|
{
|
|
// We'll get a zero cell unstructured grid which matches the input grid
|
|
tempCellList = vtkIdList::New();
|
|
}
|
|
else
|
|
{
|
|
tempCellList = cells;
|
|
}
|
|
|
|
vtkUnstructuredGrid *subGrid = vtkDistributedDataFilter::ExtractCells(
|
|
&tempCellList, 1, deleteCellLists, in, mmd);
|
|
|
|
if (tempCellList != cells)
|
|
{
|
|
tempCellList->Delete();
|
|
}
|
|
|
|
return subGrid;
|
|
}
|
|
vtkUnstructuredGrid
|
|
*vtkDistributedDataFilter::ExtractCells(vtkIdList **cells, int nlists,
|
|
int deleteCellLists, vtkDataSet *in, vtkModelMetadata *mmd)
|
|
{
|
|
vtkDataSet* tmpInput = in->NewInstance();
|
|
tmpInput->ShallowCopy(in);
|
|
|
|
vtkExtractCells *extCells = vtkExtractCells::New();
|
|
|
|
extCells->SetInput(tmpInput);
|
|
|
|
for (int i=0; i<nlists; i++)
|
|
{
|
|
if (cells[i])
|
|
{
|
|
extCells->AddCellList(cells[i]);
|
|
|
|
if (deleteCellLists)
|
|
{
|
|
cells[i]->Delete();
|
|
}
|
|
}
|
|
}
|
|
|
|
extCells->Update();
|
|
|
|
// If this process has no cells for these regions, a ugrid gets
|
|
// created anyway with field array information
|
|
|
|
vtkUnstructuredGrid *keepGrid = vtkUnstructuredGrid::New();
|
|
keepGrid->ShallowCopy(extCells->GetOutput());
|
|
|
|
extCells->Delete();
|
|
|
|
tmpInput->Delete();
|
|
|
|
if (mmd)
|
|
{
|
|
this->AddMetadata(keepGrid, mmd);
|
|
}
|
|
|
|
return keepGrid;
|
|
}
|
|
vtkUnstructuredGrid
|
|
*vtkDistributedDataFilter::ExtractZeroCellGrid(vtkDataSet *in,
|
|
vtkModelMetadata *mmd)
|
|
{
|
|
vtkDataSet* tmpInput = in->NewInstance();
|
|
tmpInput->ShallowCopy(in);
|
|
|
|
vtkExtractCells *extCells = vtkExtractCells::New();
|
|
|
|
extCells->SetInput(tmpInput);
|
|
|
|
extCells->Update(); // extract no cells
|
|
|
|
vtkUnstructuredGrid *keepGrid = vtkUnstructuredGrid::New();
|
|
keepGrid->ShallowCopy(extCells->GetOutput());
|
|
|
|
extCells->Delete();
|
|
|
|
tmpInput->Delete();
|
|
|
|
if (mmd)
|
|
{
|
|
this->AddMetadata(keepGrid, mmd);
|
|
}
|
|
|
|
return keepGrid;
|
|
}
|
|
|
|
// To save on storage, we return actual pointers into the vtkKdTree's lists
|
|
// of cell IDs. So don't free the memory they are pointing to.
|
|
// vtkKdTree::DeleteCellLists will delete them all when we're done.
|
|
|
|
vtkIdList **vtkDistributedDataFilter::GetCellIdsForProcess(int proc, int *nlists)
|
|
{
|
|
*nlists = 0;
|
|
|
|
vtkIntArray *regions = vtkIntArray::New();
|
|
|
|
int nregions = this->Kdtree->GetRegionAssignmentList(proc, regions);
|
|
|
|
if (nregions == 0)
|
|
{
|
|
return NULL;
|
|
}
|
|
|
|
*nlists = nregions;
|
|
|
|
if (this->IncludeAllIntersectingCells)
|
|
{
|
|
*nlists *= 2;
|
|
}
|
|
|
|
vtkIdList **lists = new vtkIdList * [*nlists];
|
|
|
|
int nextList = 0;
|
|
|
|
for (int reg=0; reg < nregions; reg++)
|
|
{
|
|
lists[nextList++] = this->Kdtree->GetCellList(regions->GetValue(reg));
|
|
|
|
if (this->IncludeAllIntersectingCells)
|
|
{
|
|
lists[nextList++] = this->Kdtree->GetBoundaryCellList(regions->GetValue(reg));
|
|
}
|
|
}
|
|
|
|
regions->Delete();
|
|
|
|
return lists;
|
|
}
|
|
|
|
//-------------------------------------------------------------------------
|
|
// Code related to clipping cells to the spatial region
|
|
//-------------------------------------------------------------------------
|
|
|
|
static int insideBoxFunction(vtkIdType cellId, vtkUnstructuredGrid *grid, void *data)
|
|
{
|
|
char *arrayName = (char *)data;
|
|
|
|
vtkDataArray *da= grid->GetCellData()->GetArray(arrayName);
|
|
vtkUnsignedCharArray *inside = vtkUnsignedCharArray::SafeDownCast(da);
|
|
|
|
unsigned char where = inside->GetValue(cellId);
|
|
|
|
return where; // 1 if cell is inside spatial region, 0 otherwise
|
|
}
|
|
void vtkDistributedDataFilter::AddConstantUnsignedCharPointArray(
|
|
vtkUnstructuredGrid *grid, const char *arrayName, unsigned char val)
|
|
{
|
|
vtkIdType npoints = grid->GetNumberOfPoints();
|
|
|
|
unsigned char *vals = new unsigned char [npoints];
|
|
|
|
memset(vals, val, npoints);
|
|
|
|
vtkUnsignedCharArray *Array = vtkUnsignedCharArray::New();
|
|
Array->SetName(arrayName);
|
|
Array->SetArray(vals, npoints, 0);
|
|
|
|
grid->GetPointData()->AddArray(Array);
|
|
|
|
Array->Delete();
|
|
}
|
|
void vtkDistributedDataFilter::AddConstantUnsignedCharCellArray(
|
|
vtkUnstructuredGrid *grid, const char *arrayName, unsigned char val)
|
|
{
|
|
vtkIdType ncells = grid->GetNumberOfCells();
|
|
|
|
unsigned char *vals = new unsigned char [ncells];
|
|
|
|
memset(vals, val, ncells);
|
|
|
|
vtkUnsignedCharArray *Array = vtkUnsignedCharArray::New();
|
|
Array->SetName(arrayName);
|
|
Array->SetArray(vals, ncells, 0);
|
|
|
|
grid->GetCellData()->AddArray(Array);
|
|
|
|
Array->Delete();
|
|
}
|
|
|
|
// this is here temporarily, until vtkBoxClipDataSet is fixed to
|
|
// be able to generate the clipped output
|
|
|
|
void vtkDistributedDataFilter::ClipWithVtkClipDataSet(
|
|
vtkUnstructuredGrid *grid, double *bounds,
|
|
vtkUnstructuredGrid **outside, vtkUnstructuredGrid **inside)
|
|
{
|
|
vtkUnstructuredGrid *in;
|
|
vtkUnstructuredGrid *out ;
|
|
|
|
vtkClipDataSet *clipped = vtkClipDataSet::New();
|
|
|
|
vtkBox *box = vtkBox::New();
|
|
box->SetBounds(bounds);
|
|
|
|
clipped->SetClipFunction(box);
|
|
box->Delete();
|
|
clipped->SetValue(0.0);
|
|
clipped->InsideOutOn();
|
|
|
|
clipped->SetInput(grid);
|
|
|
|
if (outside)
|
|
{
|
|
clipped->GenerateClippedOutputOn();
|
|
}
|
|
|
|
clipped->Update();
|
|
|
|
if (outside)
|
|
{
|
|
out = clipped->GetClippedOutput();
|
|
out->Register(this);
|
|
*outside = out;
|
|
}
|
|
|
|
in = clipped->GetOutput();
|
|
in->Register(this);
|
|
*inside = in;
|
|
|
|
|
|
clipped->Delete();
|
|
}
|
|
|
|
// In general, vtkBoxClipDataSet is much faster and makes fewer errors.
|
|
|
|
void vtkDistributedDataFilter::ClipWithBoxClipDataSet(
|
|
vtkUnstructuredGrid *grid, double *bounds,
|
|
vtkUnstructuredGrid **outside, vtkUnstructuredGrid **inside)
|
|
{
|
|
vtkUnstructuredGrid *in;
|
|
vtkUnstructuredGrid *out ;
|
|
|
|
vtkBoxClipDataSet *clipped = vtkBoxClipDataSet::New();
|
|
|
|
clipped->SetBoxClip(bounds[0], bounds[1],
|
|
bounds[2], bounds[3], bounds[4], bounds[5]);
|
|
|
|
clipped->SetInput(grid);
|
|
|
|
if (outside)
|
|
{
|
|
clipped->GenerateClippedOutputOn();
|
|
}
|
|
|
|
clipped->Update();
|
|
|
|
if (outside)
|
|
{
|
|
out = clipped->GetClippedOutput();
|
|
out->Register(this);
|
|
*outside = out;
|
|
}
|
|
|
|
in = clipped->GetOutput();
|
|
in->Register(this);
|
|
*inside = in;
|
|
|
|
clipped->Delete();
|
|
}
|
|
|
|
void vtkDistributedDataFilter::ClipCellsToSpatialRegion(vtkUnstructuredGrid *grid)
|
|
{
|
|
if (this->ConvexSubRegionBounds == NULL)
|
|
{
|
|
this->ComputeMyRegionBounds();
|
|
}
|
|
|
|
if (this->NumConvexSubRegions > 1)
|
|
{
|
|
// here we would need to divide the grid into a separate grid for
|
|
// each convex region, and then do the clipping
|
|
|
|
vtkErrorMacro(<<
|
|
"vtkDistributedDataFilter::ClipCellsToSpatialRegion - "
|
|
"assigned regions do not form a single convex region");
|
|
|
|
return ;
|
|
}
|
|
|
|
double *bounds = this->ConvexSubRegionBounds;
|
|
|
|
if (this->GhostLevel > 0)
|
|
{
|
|
// We need cells outside the clip box as well.
|
|
|
|
vtkUnstructuredGrid *outside;
|
|
vtkUnstructuredGrid *inside;
|
|
|
|
#if 1
|
|
this->ClipWithBoxClipDataSet(grid, bounds, &outside, &inside);
|
|
#else
|
|
this->ClipWithVtkClipDataSet(grid, bounds, &outside, &inside);
|
|
#endif
|
|
|
|
grid->Initialize();
|
|
|
|
// Mark the outside cells with a 0, the inside cells with a 1.
|
|
|
|
int arrayNameLen = strlen(TEMP_INSIDE_BOX_FLAG);
|
|
char *arrayName = new char [arrayNameLen + 1];
|
|
strcpy(arrayName, TEMP_INSIDE_BOX_FLAG);
|
|
vtkDistributedDataFilter::AddConstantUnsignedCharCellArray(outside, arrayName, 0);
|
|
vtkDistributedDataFilter::AddConstantUnsignedCharCellArray(inside, arrayName, 1);
|
|
|
|
// Combine inside and outside into a single ugrid.
|
|
|
|
vtkDataSet *grids[2];
|
|
grids[0] = inside;
|
|
grids[1] = outside;
|
|
|
|
vtkUnstructuredGrid *combined =
|
|
vtkDistributedDataFilter::MergeGrids(grids, 2, DeleteYes, NULL,
|
|
(float)this->Kdtree->GetFudgeFactor(), NULL);
|
|
|
|
// Extract the piece inside the box (level 0) and the requested
|
|
// number of levels of ghost cells.
|
|
|
|
vtkExtractUserDefinedPiece *ep = vtkExtractUserDefinedPiece::New();
|
|
|
|
ep->SetConstantData(arrayName, arrayNameLen + 1);
|
|
ep->SetPieceFunction(insideBoxFunction);
|
|
ep->CreateGhostCellsOn();
|
|
|
|
ep->GetExecutive()->GetOutputInformation(0)->Set(
|
|
vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS(),
|
|
this->GhostLevel);
|
|
ep->SetInput(combined);
|
|
|
|
ep->Update();
|
|
|
|
grid->ShallowCopy(ep->GetOutput());
|
|
grid->GetCellData()->RemoveArray(arrayName);
|
|
|
|
ep->Delete();
|
|
combined->Delete();
|
|
|
|
delete [] arrayName;
|
|
}
|
|
else
|
|
{
|
|
vtkUnstructuredGrid *inside;
|
|
|
|
#if 1
|
|
this->ClipWithBoxClipDataSet(grid, bounds, NULL, &inside);
|
|
#else
|
|
this->ClipWithVtkClipDataSet(grid, bounds, NULL, &inside);
|
|
#endif
|
|
|
|
grid->ShallowCopy(inside);
|
|
inside->Delete();
|
|
}
|
|
|
|
return;
|
|
}
|
|
|
|
//-------------------------------------------------------------------------
|
|
// Code related to assigning global node IDs and cell IDs
|
|
//-------------------------------------------------------------------------
|
|
|
|
int vtkDistributedDataFilter::AssignGlobalNodeIds(vtkUnstructuredGrid *grid)
|
|
{
|
|
int nprocs = this->NumProcesses;
|
|
int pid;
|
|
int ptId;
|
|
vtkIdType nGridPoints = grid->GetNumberOfPoints();
|
|
|
|
int *numPointsOutside = new int [nprocs];
|
|
memset(numPointsOutside, 0, sizeof(int) * nprocs);
|
|
|
|
vtkIntArray *globalIds = vtkIntArray::New();
|
|
globalIds->SetNumberOfValues(nGridPoints);
|
|
globalIds->SetName(TEMP_NODE_ID_NAME);
|
|
|
|
// 1. Count the points in grid which lie within my assigned spatial region
|
|
|
|
int myNumPointsInside = 0;
|
|
|
|
for (ptId = 0; ptId < nGridPoints; ptId++)
|
|
{
|
|
double *pt = grid->GetPoints()->GetPoint(ptId);
|
|
|
|
if (this->InMySpatialRegion(pt[0], pt[1], pt[2]))
|
|
{
|
|
globalIds->SetValue(ptId, 0); // flag it as mine
|
|
myNumPointsInside++;
|
|
}
|
|
else
|
|
{
|
|
// Well, whose region is this point in?
|
|
|
|
int regionId = this->Kdtree->GetRegionContainingPoint(pt[0],pt[1],pt[2]);
|
|
|
|
pid = this->Kdtree->GetProcessAssignedToRegion(regionId);
|
|
|
|
numPointsOutside[pid]++;
|
|
|
|
pid += 1;
|
|
pid *= -1;
|
|
|
|
globalIds->SetValue(ptId, pid); // a flag
|
|
}
|
|
}
|
|
|
|
// 2. Gather and Broadcast this number of "Inside" points for each process.
|
|
|
|
vtkIntArray *numPointsInside = this->ExchangeCounts(myNumPointsInside, 0x0013);
|
|
|
|
// 3. Assign global Ids to the points inside my spatial region
|
|
|
|
int firstId = 0;
|
|
int numGlobalIdsSoFar = 0;
|
|
|
|
for (pid = 0; pid < nprocs; pid++)
|
|
{
|
|
if (pid < this->MyId)
|
|
{
|
|
firstId += numPointsInside->GetValue(pid);
|
|
}
|
|
numGlobalIdsSoFar += numPointsInside->GetValue(pid);
|
|
}
|
|
|
|
numPointsInside->Delete();
|
|
|
|
for (ptId = 0; ptId < nGridPoints; ptId++)
|
|
{
|
|
if (globalIds->GetValue(ptId) == 0)
|
|
{
|
|
globalIds->SetValue(ptId, firstId++);
|
|
}
|
|
}
|
|
|
|
// -----------------------------------------------------------------
|
|
// All processes have assigned global IDs to the points in their grid
|
|
// which lie within their assigned spatial region.
|
|
// Now they have to get the IDs for the
|
|
// points in their grid which lie outside their region, and which
|
|
// are within the spatial region of another process.
|
|
// -----------------------------------------------------------------
|
|
|
|
// 4. For every other process, build a list of points I have
|
|
// which are in the region of that process. In practice, the
|
|
// processes for which I need to request points IDs should be
|
|
// a small subset of all the other processes.
|
|
|
|
// question: if the vtkPointArray has type double, should we
|
|
// send doubles instead of floats to insure we get the right
|
|
// global ID back?
|
|
|
|
vtkFloatArray **ptarrayOut = new vtkFloatArray * [nprocs];
|
|
memset(ptarrayOut, 0, sizeof(vtkFloatArray *) * nprocs);
|
|
|
|
vtkIntArray **localIds = new vtkIntArray * [nprocs];
|
|
memset(localIds, 0, sizeof(vtkIntArray *) * nprocs);
|
|
|
|
int *next = new int [nprocs];
|
|
int *next3 = new int [nprocs];
|
|
|
|
for (ptId = 0; ptId < nGridPoints; ptId++)
|
|
{
|
|
pid = globalIds->GetValue(ptId);
|
|
|
|
if (pid >= 0)
|
|
{
|
|
continue; // that's one of mine
|
|
}
|
|
|
|
pid *= -1;
|
|
pid -= 1;
|
|
|
|
if (ptarrayOut[pid] == NULL)
|
|
{
|
|
int npoints = numPointsOutside[pid];
|
|
|
|
ptarrayOut[pid] = vtkFloatArray::New();
|
|
ptarrayOut[pid]->SetNumberOfValues(npoints * 3);
|
|
|
|
localIds[pid] = vtkIntArray::New();
|
|
localIds[pid]->SetNumberOfValues(npoints);
|
|
|
|
next[pid] = 0;
|
|
next3[pid] = 0;
|
|
}
|
|
|
|
localIds[pid]->SetValue(next[pid]++, ptId);
|
|
|
|
double *dp = grid->GetPoints()->GetPoint(ptId);
|
|
|
|
ptarrayOut[pid]->SetValue(next3[pid]++, (float)dp[0]);
|
|
ptarrayOut[pid]->SetValue(next3[pid]++, (float)dp[1]);
|
|
ptarrayOut[pid]->SetValue(next3[pid]++, (float)dp[2]);
|
|
}
|
|
|
|
delete [] numPointsOutside;
|
|
delete [] next;
|
|
delete [] next3;
|
|
|
|
// 5. Do pairwise exchanges of the points we want global IDs for,
|
|
// and delete outgoing point arrays.
|
|
|
|
vtkFloatArray **ptarrayIn = this->ExchangeFloatArrays(ptarrayOut,
|
|
DeleteYes, 0x0014);
|
|
|
|
// 6. Find the global point IDs that have been requested of me,
|
|
// and delete incoming point arrays. Count "missing points":
|
|
// the number of unique points I receive which are not in my
|
|
// grid (this may happen if IncludeAllIntersectingCells is OFF).
|
|
|
|
int myNumMissingPoints = 0;
|
|
|
|
vtkIntArray **idarrayOut =
|
|
this->FindGlobalPointIds(ptarrayIn, globalIds, grid, myNumMissingPoints);
|
|
|
|
vtkIntArray *missingCount = this->ExchangeCounts(myNumMissingPoints, 0x0015);
|
|
|
|
if (this->IncludeAllIntersectingCells == 1)
|
|
{
|
|
// Make sure all points were found
|
|
|
|
int aok = 1;
|
|
for (pid=0; pid<nprocs; pid++)
|
|
{
|
|
if (missingCount->GetValue(pid) > 0)
|
|
{
|
|
vtkErrorMacro(<<
|
|
"vtkDistributedDataFilter::AssignGlobalNodeIds bad point");
|
|
aok = 0;
|
|
break;
|
|
}
|
|
}
|
|
if (!aok)
|
|
{
|
|
this->FreeIntArrays(idarrayOut);
|
|
this->FreeIntArrays(localIds);
|
|
missingCount->Delete();
|
|
globalIds->Delete();
|
|
|
|
return 1;
|
|
}
|
|
}
|
|
|
|
// 7. Do pairwise exchanges of the global point IDs, and delete the
|
|
// outgoing point ID arrays.
|
|
|
|
vtkIntArray **idarrayIn = this->ExchangeIntArrays(idarrayOut,
|
|
DeleteYes, 0x0016);
|
|
|
|
// 8. It's possible (if IncludeAllIntersectingCells is OFF) that some
|
|
// processes had "missing points". Process A has a point P in it's
|
|
// grid which lies in the spatial region of process B. But P is not
|
|
// in process B's grid. We need to assign global IDs to these points
|
|
// too.
|
|
|
|
int *missingId = new int [nprocs];
|
|
|
|
if (this->IncludeAllIntersectingCells == 0)
|
|
{
|
|
missingId[0] = numGlobalIdsSoFar;
|
|
|
|
for (pid = 1; pid < nprocs; pid++)
|
|
{
|
|
int prev = pid - 1;
|
|
missingId[pid] = missingId[prev] + missingCount->GetValue(prev);
|
|
}
|
|
}
|
|
|
|
missingCount->Delete();
|
|
|
|
// 9. Update my ugrid with these mutually agreed upon global point IDs
|
|
|
|
for (pid = 0; pid < nprocs; pid++)
|
|
{
|
|
if (idarrayIn[pid] == NULL)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
int count = idarrayIn[pid]->GetNumberOfTuples();
|
|
|
|
for (ptId = 0; ptId < count; ptId++)
|
|
{
|
|
int myLocalId = localIds[pid]->GetValue(ptId);
|
|
int yourGlobalId = idarrayIn[pid]->GetValue(ptId);
|
|
|
|
if (yourGlobalId >= 0)
|
|
{
|
|
globalIds->SetValue(myLocalId, yourGlobalId);
|
|
}
|
|
else
|
|
{
|
|
int ptIdOffset = yourGlobalId * -1;
|
|
ptIdOffset -= 1;
|
|
|
|
globalIds->SetValue(myLocalId, missingId[pid] + ptIdOffset);
|
|
}
|
|
}
|
|
localIds[pid]->Delete();
|
|
idarrayIn[pid]->Delete();
|
|
}
|
|
|
|
delete [] localIds;
|
|
delete [] idarrayIn;
|
|
delete [] missingId;
|
|
|
|
grid->GetPointData()->AddArray(globalIds);
|
|
globalIds->Delete();
|
|
|
|
this->SetGlobalNodeIdArrayName(TEMP_NODE_ID_NAME);
|
|
|
|
return 0;
|
|
}
|
|
|
|
// If grids were distributed with IncludeAllIntersectingCells OFF, it's
|
|
// possible there are points in my spatial region that are not in my
|
|
// grid. They need global Ids, so I will keep track of how many such unique
|
|
// points I receive from other processes, and will assign them temporary
|
|
// IDs. They will get permanent IDs later on.
|
|
|
|
vtkIntArray **vtkDistributedDataFilter::FindGlobalPointIds(
|
|
vtkFloatArray **ptarray, vtkIntArray *ids, vtkUnstructuredGrid *grid,
|
|
int &numUniqueMissingPoints)
|
|
{
|
|
int nprocs = this->NumProcesses;
|
|
vtkIntArray **gids = new vtkIntArray * [nprocs];
|
|
|
|
if (grid->GetNumberOfCells() == 0)
|
|
{
|
|
// There are no cells in my assigned region
|
|
|
|
memset(gids, 0, sizeof(vtkIntArray *) * nprocs);
|
|
|
|
return gids;
|
|
}
|
|
|
|
vtkKdTree *kd = vtkKdTree::New();
|
|
|
|
kd->BuildLocatorFromPoints(grid->GetPoints());
|
|
|
|
int procId;
|
|
int ptId, localId;
|
|
|
|
vtkPointLocator *pl = NULL;
|
|
vtkPoints *missingPoints = NULL;
|
|
|
|
if (this->IncludeAllIntersectingCells == 0)
|
|
{
|
|
if (this->ConvexSubRegionBounds == NULL)
|
|
{
|
|
this->ComputeMyRegionBounds();
|
|
}
|
|
pl = vtkPointLocator::New();
|
|
pl->SetTolerance(this->Kdtree->GetFudgeFactor());
|
|
missingPoints = vtkPoints::New();
|
|
pl->InitPointInsertion(missingPoints, this->ConvexSubRegionBounds);
|
|
}
|
|
|
|
for (procId = 0; procId < nprocs; procId++)
|
|
{
|
|
if ((ptarray[procId] == NULL) ||
|
|
(ptarray[procId]->GetNumberOfTuples() == 0))
|
|
{
|
|
gids[procId] = NULL;
|
|
if (ptarray[procId]) ptarray[procId]->Delete();
|
|
continue;
|
|
}
|
|
|
|
gids[procId] = vtkIntArray::New();
|
|
|
|
int npoints = ptarray[procId]->GetNumberOfTuples() / 3;
|
|
|
|
gids[procId]->SetNumberOfValues(npoints);
|
|
int next = 0;
|
|
|
|
float *pt = ptarray[procId]->GetPointer(0);
|
|
|
|
for (ptId = 0; ptId < npoints; ptId++)
|
|
{
|
|
localId = kd->FindPoint((double)pt[0], (double)pt[1], (double)pt[2]);
|
|
|
|
if (localId >= 0)
|
|
{
|
|
gids[procId]->SetValue(next++, ids->GetValue(localId)); // global Id
|
|
}
|
|
else
|
|
{
|
|
// This point is not in my grid
|
|
|
|
if (this->IncludeAllIntersectingCells)
|
|
{
|
|
// This is an error
|
|
gids[procId]->SetValue(next++, -1);
|
|
numUniqueMissingPoints++;
|
|
}
|
|
else
|
|
{
|
|
// Flag these with a negative point ID. We'll assign
|
|
// them real point IDs later.
|
|
|
|
vtkIdType nextId;
|
|
double dpt[3];
|
|
dpt[0] = pt[0]; dpt[1] = pt[1]; dpt[2] = pt[2];
|
|
pl->InsertUniquePoint(dpt, nextId);
|
|
|
|
nextId += 1;
|
|
nextId *= -1;
|
|
gids[procId]->SetValue(next++, nextId);
|
|
}
|
|
}
|
|
pt += 3;
|
|
}
|
|
|
|
ptarray[procId]->Delete();
|
|
}
|
|
|
|
delete [] ptarray;
|
|
|
|
kd->Delete();
|
|
|
|
if (missingPoints)
|
|
{
|
|
numUniqueMissingPoints = missingPoints->GetNumberOfPoints();
|
|
missingPoints->Delete();
|
|
pl->Delete();
|
|
}
|
|
|
|
return gids;
|
|
}
|
|
int vtkDistributedDataFilter::AssignGlobalElementIds(vtkDataSet *in)
|
|
{
|
|
int i;
|
|
int myNumCells = in->GetNumberOfCells();
|
|
vtkIntArray *numCells = this->ExchangeCounts(myNumCells, 0x0017);
|
|
|
|
vtkIntArray *globalCellIds = vtkIntArray::New();
|
|
globalCellIds->SetNumberOfValues(myNumCells);
|
|
globalCellIds->SetName(TEMP_ELEMENT_ID_NAME);
|
|
|
|
int StartId = 0;
|
|
|
|
for (i=0; i < this->MyId; i++)
|
|
{
|
|
StartId += numCells->GetValue(i);
|
|
}
|
|
|
|
numCells->Delete();
|
|
|
|
for (i=0; i<myNumCells; i++)
|
|
{
|
|
globalCellIds->SetValue(i, StartId++);
|
|
}
|
|
|
|
in->GetCellData()->AddArray(globalCellIds);
|
|
|
|
globalCellIds->Delete();
|
|
|
|
this->SetGlobalElementIdArrayName(TEMP_ELEMENT_ID_NAME);
|
|
|
|
return 0;
|
|
}
|
|
|
|
//-------------------------------------------------------------------------
|
|
// Code related to acquiring ghost cells
|
|
//-------------------------------------------------------------------------
|
|
|
|
int vtkDistributedDataFilter::InMySpatialRegion(float x, float y, float z)
|
|
{
|
|
return this->InMySpatialRegion((double)x, (double)y, (double)z);
|
|
}
|
|
int vtkDistributedDataFilter::InMySpatialRegion(double x, double y, double z)
|
|
{
|
|
if (this->ConvexSubRegionBounds == NULL)
|
|
{
|
|
this->ComputeMyRegionBounds();
|
|
}
|
|
|
|
double *box = this->ConvexSubRegionBounds;
|
|
|
|
if (!box)
|
|
{
|
|
return 0;
|
|
}
|
|
|
|
// To avoid ambiguity, a point on a boundary is assigned to
|
|
// the region for which it is on the upper boundary. Or
|
|
// (in one dimension) the region between points A and B
|
|
// contains all points p such that A < p <= B.
|
|
|
|
if ( (x <= box[0]) || (x > box[1]) ||
|
|
(y <= box[2]) || (y > box[3]) ||
|
|
(z <= box[4]) || (z > box[5]) )
|
|
{
|
|
return 0;
|
|
}
|
|
|
|
return 1;
|
|
}
|
|
int vtkDistributedDataFilter::StrictlyInsideMyBounds(float x, float y, float z)
|
|
{
|
|
return this->StrictlyInsideMyBounds((double)x, (double)y, (double)z);
|
|
}
|
|
int vtkDistributedDataFilter::StrictlyInsideMyBounds(double x, double y, double z)
|
|
{
|
|
if (this->ConvexSubRegionBounds == NULL)
|
|
{
|
|
this->ComputeMyRegionBounds();
|
|
}
|
|
|
|
double *box = this->ConvexSubRegionBounds;
|
|
|
|
if (!box)
|
|
{
|
|
return 0;
|
|
}
|
|
|
|
if ( (x <= box[0]) || (x >= box[1]) ||
|
|
(y <= box[2]) || (y >= box[3]) ||
|
|
(z <= box[4]) || (z >= box[5]) )
|
|
{
|
|
return 0;
|
|
}
|
|
|
|
return 1;
|
|
}
|
|
|
|
vtkIntArray **vtkDistributedDataFilter::MakeProcessLists(
|
|
vtkIntArray **pointIds,
|
|
vtkDistributedDataFilterSTLCloak *procs)
|
|
{
|
|
// Build a list of pointId/processId pairs for each process that
|
|
// sent me point IDs. The process Ids are all those processes
|
|
// that had the specified point in their ghost level zero grid.
|
|
|
|
int nprocs = this->NumProcesses;
|
|
|
|
vtkstd::multimap<int, int>::iterator mapIt;
|
|
|
|
vtkIntArray **processList = new vtkIntArray * [nprocs];
|
|
memset(processList, 0, sizeof (vtkIntArray *) * nprocs);
|
|
|
|
for (int i=0; i<nprocs; i++)
|
|
{
|
|
if (pointIds[i] == NULL)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
int size = pointIds[i]->GetNumberOfTuples();
|
|
|
|
if (size > 0)
|
|
{
|
|
for (int j=0; j<size; )
|
|
{
|
|
// These are all the points in my spatial region
|
|
// for which process "i" needs ghost cells.
|
|
|
|
int gid = pointIds[i]->GetValue(j);
|
|
int ncells = pointIds[i]->GetValue(j+1);
|
|
|
|
mapIt = procs->IntMultiMap.find(gid);
|
|
|
|
if (mapIt != procs->IntMultiMap.end())
|
|
{
|
|
while (mapIt->first == gid)
|
|
{
|
|
int processId = mapIt->second;
|
|
|
|
if (processId != i)
|
|
{
|
|
// Process "i" needs to know that process
|
|
// "processId" also has cells using this point
|
|
|
|
if (processList[i] == NULL)
|
|
{
|
|
processList[i] = vtkIntArray::New();
|
|
}
|
|
processList[i]->InsertNextValue(gid);
|
|
processList[i]->InsertNextValue(processId);
|
|
}
|
|
++mapIt;
|
|
}
|
|
}
|
|
j += (2 + ncells);
|
|
}
|
|
}
|
|
}
|
|
|
|
return processList;
|
|
}
|
|
vtkIntArray *vtkDistributedDataFilter::AddPointAndCells(
|
|
int gid, int localId, vtkUnstructuredGrid *grid,
|
|
int *gidCells, vtkIntArray *ids)
|
|
{
|
|
if (ids == NULL)
|
|
{
|
|
ids = vtkIntArray::New();
|
|
}
|
|
|
|
ids->InsertNextValue(gid);
|
|
|
|
vtkIdList *cellList = vtkIdList::New();
|
|
|
|
grid->GetPointCells(localId, cellList);
|
|
|
|
vtkIdType numCells = cellList->GetNumberOfIds();
|
|
|
|
ids->InsertNextValue((int)numCells);
|
|
|
|
for (int j=0; j<numCells; j++)
|
|
{
|
|
int globalCellId = gidCells[cellList->GetId(j)];
|
|
ids->InsertNextValue(globalCellId);
|
|
}
|
|
|
|
cellList->Delete();
|
|
|
|
return ids;
|
|
}
|
|
vtkIntArray **vtkDistributedDataFilter::GetGhostPointIds(
|
|
int ghostLevel, vtkUnstructuredGrid *grid,
|
|
int AddCellsIAlreadyHave)
|
|
{
|
|
int nprocs = this->NumProcesses;
|
|
int me = this->MyId;
|
|
vtkIdType numPoints = grid->GetNumberOfPoints();
|
|
|
|
vtkIntArray **ghostPtIds = new vtkIntArray * [nprocs];
|
|
memset(ghostPtIds, 0, sizeof(vtkIntArray *) * nprocs);
|
|
|
|
if (numPoints < 1)
|
|
{
|
|
return ghostPtIds;
|
|
}
|
|
|
|
int processId = -1;
|
|
int regionId = -1;
|
|
|
|
vtkPKdTree *kd = this->Kdtree;
|
|
|
|
vtkPoints *pts = grid->GetPoints();
|
|
|
|
int *gidsPoint = this->GetGlobalNodeIds(grid);
|
|
int *gidsCell = this->GetGlobalElementIds(grid);
|
|
|
|
|
|
vtkDataArray *da = grid->GetPointData()->GetArray("vtkGhostLevels");
|
|
vtkUnsignedCharArray *uca = vtkUnsignedCharArray::SafeDownCast(da);
|
|
unsigned char *levels = uca->GetPointer(0);
|
|
|
|
unsigned char level = (unsigned char)(ghostLevel - 1);
|
|
|
|
for (int i=0; i<numPoints; i++)
|
|
{
|
|
double *pt = pts->GetPoint(i);
|
|
regionId = kd->GetRegionContainingPoint(pt[0], pt[1], pt[2]);
|
|
processId = kd->GetProcessAssignedToRegion(regionId);
|
|
|
|
if (ghostLevel == 1)
|
|
{
|
|
// I want all points that are outside my spatial region
|
|
|
|
if (processId == me)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
// Don't include points that are not part of any cell
|
|
|
|
int used = vtkDistributedDataFilter::LocalPointIdIsUsed(grid, i);
|
|
|
|
if (!used)
|
|
{
|
|
continue;
|
|
}
|
|
}
|
|
else if (levels[i] != level)
|
|
{
|
|
continue; // I want all points having the correct ghost level
|
|
}
|
|
|
|
int gid = gidsPoint[i];
|
|
|
|
if (AddCellsIAlreadyHave)
|
|
{
|
|
// To speed up exchange of ghost cells and creation of
|
|
// new ghost cell grid, we tell other
|
|
// processes which cells we already have, so they don't
|
|
// send them to us.
|
|
|
|
ghostPtIds[processId] =
|
|
vtkDistributedDataFilter::AddPointAndCells(gid, i, grid, gidsCell,
|
|
ghostPtIds[processId]);
|
|
}
|
|
else
|
|
{
|
|
if (ghostPtIds[processId] == NULL)
|
|
{
|
|
ghostPtIds[processId] = vtkIntArray::New();
|
|
}
|
|
ghostPtIds[processId]->InsertNextValue(gid);
|
|
ghostPtIds[processId]->InsertNextValue(0);
|
|
}
|
|
}
|
|
return ghostPtIds;
|
|
}
|
|
int vtkDistributedDataFilter::LocalPointIdIsUsed(
|
|
vtkUnstructuredGrid *grid, int ptId)
|
|
{
|
|
int used = 1;
|
|
|
|
int numPoints = grid->GetNumberOfPoints();
|
|
|
|
if ((ptId < 0) || (ptId >= numPoints))
|
|
{
|
|
used = 0;
|
|
}
|
|
else
|
|
{
|
|
vtkIdType id = (vtkIdType)ptId;
|
|
vtkIdList *cellList = vtkIdList::New();
|
|
|
|
grid->GetPointCells(id, cellList);
|
|
|
|
if (cellList->GetNumberOfIds() == 0)
|
|
{
|
|
used = 0;
|
|
}
|
|
|
|
cellList->Delete();
|
|
}
|
|
|
|
return used;
|
|
}
|
|
int vtkDistributedDataFilter::GlobalPointIdIsUsed(vtkUnstructuredGrid *grid,
|
|
int ptId, vtkDistributedDataFilterSTLCloak *globalToLocal)
|
|
{
|
|
int used = 1;
|
|
|
|
vtkstd::map<int, int>::iterator mapIt;
|
|
|
|
mapIt = globalToLocal->IntMap.find(ptId);
|
|
|
|
if (mapIt == globalToLocal->IntMap.end())
|
|
{
|
|
used = 0;
|
|
}
|
|
else
|
|
{
|
|
int id = mapIt->second;
|
|
|
|
used = vtkDistributedDataFilter::LocalPointIdIsUsed(grid, id);
|
|
}
|
|
|
|
return used;
|
|
}
|
|
int vtkDistributedDataFilter::FindId(vtkIntArray *ids, int gid, int startLoc)
|
|
{
|
|
int gidLoc = -1;
|
|
|
|
if (ids == NULL)
|
|
{
|
|
return gidLoc;
|
|
}
|
|
|
|
int numIds = ids->GetNumberOfTuples();
|
|
|
|
while ((ids->GetValue(startLoc) != gid) && (startLoc < numIds))
|
|
{
|
|
int ncells = ids->GetValue(++startLoc);
|
|
startLoc += (ncells + 1);
|
|
}
|
|
|
|
if (startLoc < numIds)
|
|
{
|
|
gidLoc = startLoc;
|
|
}
|
|
|
|
return gidLoc;
|
|
}
|
|
|
|
// We create an expanded grid with the required number of ghost
|
|
// cells. This is for the case where IncludeAllIntersectingCells is OFF.
|
|
// This means that when the grid was redistributed, each cell was
|
|
// uniquely assigned to one process, the process owning the spatial
|
|
// region that the cell's centroid lies in.
|
|
|
|
vtkUnstructuredGrid *
|
|
vtkDistributedDataFilter::AddGhostCellsUniqueCellAssignment(
|
|
vtkUnstructuredGrid *myGrid,
|
|
vtkDistributedDataFilterSTLCloak *globalToLocalMap)
|
|
{
|
|
int i,j,k;
|
|
int ncells=0;
|
|
int processId=0;
|
|
int gid=0;
|
|
int size=0;
|
|
|
|
int nprocs = this->NumProcesses;
|
|
int me = this->MyId;
|
|
|
|
int gl = 1;
|
|
|
|
// For each ghost level, processes request and send ghost cells
|
|
|
|
vtkUnstructuredGrid *newGhostCellGrid = NULL;
|
|
vtkIntArray **ghostPointIds = NULL;
|
|
|
|
vtkDistributedDataFilterSTLCloak *insidePointMap =
|
|
new vtkDistributedDataFilterSTLCloak;
|
|
vtkstd::multimap<int, int>::iterator mapIt;
|
|
|
|
while (gl <= this->GhostLevel)
|
|
{
|
|
// For ghost level 1, create a list for each process (not
|
|
// including me) of all points I have in that process'
|
|
// assigned region. We use this list for two purposes:
|
|
// (1) to build a list on each process of all other processes
|
|
// that have cells containing points in our region, (2)
|
|
// these are some of the points that we need ghost cells for.
|
|
//
|
|
// For ghost level above 1, create a list for each process
|
|
// (including me) of all my points in that process' assigned
|
|
// region for which I need ghost cells.
|
|
|
|
if (gl == 1)
|
|
{
|
|
ghostPointIds = this->GetGhostPointIds(gl, myGrid, 0);
|
|
}
|
|
else
|
|
{
|
|
ghostPointIds = this->GetGhostPointIds(gl, newGhostCellGrid, 1);
|
|
}
|
|
|
|
// Exchange these lists.
|
|
|
|
vtkIntArray **insideIds =
|
|
this->ExchangeIntArrays(ghostPointIds, DeleteNo,
|
|
0x0018);
|
|
|
|
if (gl == 1)
|
|
{
|
|
// For every point in my region that was sent to me by another process,
|
|
// I now know the identity of all processes having cells containing
|
|
// that point. Begin by building a mapping from point IDs to the IDs
|
|
// of processes that sent me that point.
|
|
|
|
for (i=0; i<nprocs; i++)
|
|
{
|
|
if (insideIds[i] == NULL)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
size = insideIds[i]->GetNumberOfTuples();
|
|
|
|
if (size > 0)
|
|
{
|
|
for (j=0; j<size; j+=2)
|
|
{
|
|
// map global point id to process ids
|
|
const int id = insideIds[i]->GetValue(j);
|
|
insidePointMap->IntMultiMap.insert(vtkstd::pair<const int, int>(id, i));
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Build a list of pointId/processId pairs for each process that
|
|
// sent me point IDs. To process P, for every point ID sent to me
|
|
// by P, I send the ID of every other process (not including myself
|
|
// and P) that has cells in it's ghost level 0 grid which use
|
|
// this point.
|
|
|
|
vtkIntArray **processListSent
|
|
= this->MakeProcessLists(insideIds, insidePointMap);
|
|
|
|
// Exchange these new lists.
|
|
|
|
vtkIntArray **processList =
|
|
this->ExchangeIntArrays(processListSent, DeleteYes,
|
|
0x0019);
|
|
|
|
// I now know the identity of every process having cells containing
|
|
// points I need ghost cells for. Create a request to each process
|
|
// for these cells.
|
|
|
|
vtkIntArray **ghostCellsPlease = new vtkIntArray * [nprocs];
|
|
for (i=0; i<nprocs; i++)
|
|
{
|
|
ghostCellsPlease[i] = vtkIntArray::New();
|
|
ghostCellsPlease[i]->SetNumberOfComponents(1);
|
|
}
|
|
|
|
for (i=0; i<nprocs; i++)
|
|
{
|
|
if (i == me)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (ghostPointIds[i]) // points I have in your spatial region,
|
|
{ // maybe you have cells that use them?
|
|
|
|
for (j=0; j<ghostPointIds[i]->GetNumberOfTuples(); j++)
|
|
{
|
|
ghostCellsPlease[i]->InsertNextValue(ghostPointIds[i]->GetValue(j));
|
|
}
|
|
}
|
|
if (processList[i]) // other processes you say that also have
|
|
{ // cells using those points
|
|
size = processList[i]->GetNumberOfTuples();
|
|
int *array = processList[i]->GetPointer(0);
|
|
int nextLoc = 0;
|
|
|
|
for (j=0; j < size; j += 2)
|
|
{
|
|
gid = array[j];
|
|
processId = array[j+1];
|
|
|
|
ghostCellsPlease[processId]->InsertNextValue(gid);
|
|
|
|
if (gl > 1)
|
|
{
|
|
// add the list of cells I already have for this point
|
|
|
|
int where =
|
|
vtkDistributedDataFilter::FindId(ghostPointIds[i], gid, nextLoc);
|
|
|
|
if (where < 0)
|
|
{
|
|
// error really, not sure what to do
|
|
nextLoc = 0;
|
|
ghostCellsPlease[processId]->InsertNextValue(0);
|
|
continue;
|
|
}
|
|
|
|
ncells = ghostPointIds[i]->GetValue(where + 1);
|
|
|
|
ghostCellsPlease[processId]->InsertNextValue(ncells);
|
|
|
|
for (k=0; k <ncells; k++)
|
|
{
|
|
int cellId = ghostPointIds[i]->GetValue(where + 2 + k);
|
|
ghostCellsPlease[processId]->InsertNextValue(cellId);
|
|
}
|
|
|
|
nextLoc = where;
|
|
}
|
|
else
|
|
{
|
|
ghostCellsPlease[processId]->InsertNextValue(0);
|
|
}
|
|
}
|
|
}
|
|
if ((gl==1) && insideIds[i]) // points you have in my spatial region,
|
|
{ // which I may need ghost cells for
|
|
for (j=0; j<insideIds[i]->GetNumberOfTuples();)
|
|
{
|
|
gid = insideIds[i]->GetValue(j);
|
|
int used = vtkDistributedDataFilter::GlobalPointIdIsUsed(
|
|
myGrid, gid, globalToLocalMap);
|
|
if (used)
|
|
{
|
|
ghostCellsPlease[i]->InsertNextValue(gid);
|
|
ghostCellsPlease[i]->InsertNextValue(0);
|
|
}
|
|
|
|
ncells = insideIds[i]->GetValue(j+1);
|
|
j += (ncells + 2);
|
|
}
|
|
}
|
|
}
|
|
|
|
if (gl > 1)
|
|
{
|
|
if (ghostPointIds[me]) // these points are actually inside my region
|
|
{
|
|
size = ghostPointIds[me]->GetNumberOfTuples();
|
|
|
|
for (i=0; i<size;)
|
|
{
|
|
gid = ghostPointIds[me]->GetValue(i);
|
|
ncells = ghostPointIds[me]->GetValue(i+1);
|
|
|
|
mapIt = insidePointMap->IntMultiMap.find(gid);
|
|
|
|
if (mapIt != insidePointMap->IntMultiMap.end())
|
|
{
|
|
while (mapIt->first == gid)
|
|
{
|
|
processId = mapIt->second;
|
|
ghostCellsPlease[processId]->InsertNextValue(gid);
|
|
ghostCellsPlease[processId]->InsertNextValue(ncells);
|
|
|
|
for (k=0; k<ncells; k++)
|
|
{
|
|
int cellId = ghostPointIds[me]->GetValue(i+1+k);
|
|
ghostCellsPlease[processId]->InsertNextValue(cellId);
|
|
}
|
|
|
|
++mapIt;
|
|
}
|
|
}
|
|
i += (ncells + 2);
|
|
}
|
|
}
|
|
}
|
|
|
|
this->FreeIntArrays(ghostPointIds);
|
|
this->FreeIntArrays(insideIds);
|
|
this->FreeIntArrays(processList);
|
|
|
|
// Exchange these ghost cell requests.
|
|
|
|
vtkIntArray **ghostCellRequest
|
|
= this->ExchangeIntArrays(ghostCellsPlease, DeleteYes,
|
|
0x001a);
|
|
|
|
// Build a list of cell IDs satisfying each request received.
|
|
// Delete request arrays.
|
|
|
|
vtkIdList **sendCellList =
|
|
this->BuildRequestedGrids(ghostCellRequest, myGrid, globalToLocalMap);
|
|
|
|
// Build subgrids and exchange them
|
|
|
|
vtkUnstructuredGrid *incomingGhostCells = this->ExchangeMergeSubGrids(
|
|
sendCellList, DeleteYes, myGrid, DeleteNo, DuplicateCellsNo,
|
|
GhostCellsYes, 0x001b);
|
|
|
|
delete [] sendCellList;
|
|
|
|
// Set ghost level of new cells, and merge into grid of other
|
|
// ghost cells received.
|
|
|
|
newGhostCellGrid = this->SetMergeGhostGrid(newGhostCellGrid,
|
|
incomingGhostCells, gl, globalToLocalMap);
|
|
|
|
this->UpdateProgress(this->NextProgressStep++ * this->ProgressIncrement);
|
|
|
|
gl++;
|
|
}
|
|
|
|
delete insidePointMap;
|
|
|
|
vtkUnstructuredGrid *newGrid = NULL;
|
|
|
|
if (newGhostCellGrid && (newGhostCellGrid->GetNumberOfCells() > 0))
|
|
{
|
|
vtkDataSet *grids[2];
|
|
|
|
grids[0] = myGrid;
|
|
grids[1] = newGhostCellGrid;
|
|
|
|
const char *nodeIds = this->GetGlobalNodeIdArrayName(myGrid);
|
|
|
|
newGrid =
|
|
vtkDistributedDataFilter::MergeGrids(grids, 2, DeleteYes, nodeIds, 0, NULL);
|
|
}
|
|
else
|
|
{
|
|
newGrid = myGrid;
|
|
}
|
|
|
|
return newGrid;
|
|
}
|
|
|
|
// We create an expanded grid that contains the ghost cells we need.
|
|
// This is in the case where IncludeAllIntersectingCells is ON. This
|
|
// is easier in some respects because we know if that if a point lies
|
|
// in a region owned by a particular process, that process has all
|
|
// cells which use that point. So it is easy to find ghost cells.
|
|
// On the otherhand, because cells are not uniquely assigned to regions,
|
|
// we may get multiple processes sending us the same cell, so we
|
|
// need to filter these out.
|
|
|
|
vtkUnstructuredGrid *
|
|
vtkDistributedDataFilter::AddGhostCellsDuplicateCellAssignment(
|
|
vtkUnstructuredGrid *myGrid,
|
|
vtkDistributedDataFilterSTLCloak *globalToLocalMap)
|
|
{
|
|
int i,j;
|
|
|
|
int nprocs = this->NumProcesses;
|
|
int me = this->MyId;
|
|
|
|
int gl = 1;
|
|
|
|
// For each ghost level, processes request and send ghost cells
|
|
|
|
vtkUnstructuredGrid *newGhostCellGrid = NULL;
|
|
vtkIntArray **ghostPointIds = NULL;
|
|
vtkIntArray **extraGhostPointIds = NULL;
|
|
|
|
vtkstd::map<int, int>::iterator mapIt;
|
|
|
|
vtkPoints *pts = myGrid->GetPoints();
|
|
|
|
while (gl <= this->GhostLevel)
|
|
{
|
|
// For ghost level 1, create a list for each process of points
|
|
// in my grid which lie in that other process' spatial region.
|
|
// This is normally all the points for which I need ghost cells,
|
|
// with one EXCEPTION. If a cell is axis-aligned, and a face of
|
|
// the cell is on my upper boundary, then the vertices of this
|
|
// face are in my spatial region, but I need their ghost cells.
|
|
// I can detect this case when the process across the boundary
|
|
// sends me a request for ghost cells of these points.
|
|
//
|
|
// For ghost level above 1, create a list for each process of
|
|
// points in my ghost grid which are in that process' spatial
|
|
// region and for which I need ghost cells.
|
|
|
|
if (gl == 1)
|
|
{
|
|
ghostPointIds = this->GetGhostPointIds(gl, myGrid, 1);
|
|
}
|
|
else
|
|
{
|
|
ghostPointIds = this->GetGhostPointIds(gl, newGhostCellGrid, 1);
|
|
}
|
|
|
|
// Exchange these lists.
|
|
|
|
vtkIntArray **insideIds =
|
|
this->ExchangeIntArrays(ghostPointIds, DeleteYes, 0x001c);
|
|
|
|
// For ghost level 1, examine the points Ids I received from
|
|
// other processes, to see if the exception described above
|
|
// applies and I need ghost cells from them for those points.
|
|
|
|
if (gl == 1)
|
|
{
|
|
int *gidsCell = this->GetGlobalElementIds(myGrid);
|
|
|
|
extraGhostPointIds = new vtkIntArray * [nprocs];
|
|
|
|
for (i=0; i<nprocs; i++)
|
|
{
|
|
extraGhostPointIds[i] = NULL;
|
|
|
|
if (i == me)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (insideIds[i] == NULL)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
int size = insideIds[i]->GetNumberOfTuples();
|
|
|
|
for (j=0; j<size;)
|
|
{
|
|
int gid = insideIds[i]->GetValue(j);
|
|
int ncells = insideIds[i]->GetValue(j+1);
|
|
j += (ncells + 2);
|
|
|
|
mapIt = globalToLocalMap->IntMap.find(gid);
|
|
|
|
if (mapIt == globalToLocalMap->IntMap.end())
|
|
{
|
|
// This point must be right on my boundary, and
|
|
// not connected to any cell intersecting my region.
|
|
|
|
continue;
|
|
}
|
|
int localId = mapIt->second;
|
|
|
|
double *pt = pts->GetPoint(localId);
|
|
|
|
int interior = this->StrictlyInsideMyBounds(pt[0], pt[1], pt[2]);
|
|
|
|
if (!interior)
|
|
{
|
|
extraGhostPointIds[i] = this->AddPointAndCells(gid, localId,
|
|
myGrid, gidsCell, extraGhostPointIds[i]);
|
|
}
|
|
}
|
|
}
|
|
|
|
// Exchange these lists.
|
|
|
|
vtkIntArray **extraInsideIds =
|
|
this->ExchangeIntArrays(extraGhostPointIds, DeleteYes, 0x001d);
|
|
|
|
// Add the extra point ids to the previous list
|
|
|
|
for (i=0; i<nprocs; i++)
|
|
{
|
|
if (i == me)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (extraInsideIds[i])
|
|
{
|
|
int size = extraInsideIds[i]->GetNumberOfTuples();
|
|
|
|
if (insideIds[i] == NULL)
|
|
{
|
|
insideIds[i] = vtkIntArray::New();
|
|
}
|
|
|
|
for (j=0; j<size; j++)
|
|
{
|
|
insideIds[i]->InsertNextValue(extraInsideIds[i]->GetValue(j));
|
|
}
|
|
}
|
|
}
|
|
this->FreeIntArrays(extraInsideIds);
|
|
}
|
|
|
|
// Build a list of cell IDs satisfying each request received.
|
|
|
|
vtkIdList **sendCellList =
|
|
this->BuildRequestedGrids(insideIds, myGrid, globalToLocalMap);
|
|
|
|
// Build subgrids and exchange them
|
|
|
|
vtkUnstructuredGrid *incomingGhostCells =
|
|
this->ExchangeMergeSubGrids( sendCellList, DeleteYes, myGrid, DeleteNo,
|
|
DuplicateCellsYes, GhostCellsYes, 0x001e);
|
|
|
|
delete [] sendCellList;
|
|
|
|
// Set ghost level of new cells, and merge into grid of other
|
|
// ghost cells received.
|
|
|
|
newGhostCellGrid = this->SetMergeGhostGrid(newGhostCellGrid,
|
|
incomingGhostCells, gl, globalToLocalMap);
|
|
|
|
this->UpdateProgress(this->NextProgressStep++ * this->ProgressIncrement);
|
|
|
|
gl++;
|
|
}
|
|
|
|
vtkUnstructuredGrid *newGrid = NULL;
|
|
|
|
if (newGhostCellGrid && (newGhostCellGrid->GetNumberOfCells() > 0))
|
|
{
|
|
vtkDataSet *grids[2];
|
|
|
|
grids[0] = myGrid;
|
|
grids[1] = newGhostCellGrid;
|
|
|
|
const char *nodeIds = this->GetGlobalNodeIdArrayName(myGrid);
|
|
|
|
newGrid =
|
|
vtkDistributedDataFilter::MergeGrids(grids, 2, DeleteYes, nodeIds, 0, NULL);
|
|
}
|
|
else
|
|
{
|
|
newGrid = myGrid;
|
|
}
|
|
|
|
return newGrid;
|
|
}
|
|
|
|
// For every process that sent me a list of point IDs, create a list
|
|
// of all the cells I have in my original grid containing those points.
|
|
// We omit cells the remote process already has.
|
|
|
|
vtkIdList **vtkDistributedDataFilter::BuildRequestedGrids(
|
|
vtkIntArray **globalPtIds,
|
|
vtkUnstructuredGrid *grid,
|
|
vtkDistributedDataFilterSTLCloak *ptIdMap)
|
|
{
|
|
int id, proc;
|
|
int nprocs = this->NumProcesses;
|
|
vtkIdType cellId;
|
|
vtkIdType nelts;
|
|
|
|
// for each process, create a list of the ids of cells I need
|
|
// to send to it
|
|
|
|
vtkstd::map<int, int>::iterator imap;
|
|
|
|
vtkIdList *cellList = vtkIdList::New();
|
|
|
|
vtkIdList **sendCells = new vtkIdList * [nprocs];
|
|
|
|
for (proc = 0; proc < nprocs; proc++)
|
|
{
|
|
sendCells[proc] = vtkIdList::New();
|
|
|
|
if (globalPtIds[proc] == NULL)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if ((nelts = globalPtIds[proc]->GetNumberOfTuples()) == 0)
|
|
{
|
|
globalPtIds[proc]->Delete();
|
|
continue;
|
|
}
|
|
|
|
int *ptarray = globalPtIds[proc]->GetPointer(0);
|
|
|
|
vtkstd::set<vtkIdType> subGridCellIds;
|
|
|
|
int nYourCells = 0;
|
|
|
|
for (id = 0; id < nelts; id += (nYourCells + 2))
|
|
{
|
|
int ptId = ptarray[id];
|
|
nYourCells = ptarray[id+1];
|
|
|
|
imap = ptIdMap->IntMap.find(ptId);
|
|
|
|
if (imap == ptIdMap->IntMap.end())
|
|
{
|
|
continue; // I don't have this point
|
|
}
|
|
|
|
vtkIdType myPtId = (vtkIdType)imap->second; // convert to my local point Id
|
|
|
|
grid->GetPointCells(myPtId, cellList);
|
|
|
|
vtkIdType nMyCells = cellList->GetNumberOfIds();
|
|
|
|
if (nMyCells == 0)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
if (nYourCells > 0)
|
|
{
|
|
// We don't send cells the remote process tells us it already
|
|
// has. This is much faster than removing duplicate cells on
|
|
// the receive side.
|
|
|
|
int *remoteCells = ptarray + id + 2;
|
|
int *gidCells = this->GetGlobalElementIds(grid);
|
|
|
|
vtkDistributedDataFilter::RemoveRemoteCellsFromList(cellList,
|
|
gidCells, remoteCells, nYourCells);
|
|
}
|
|
|
|
vtkIdType nSendCells = cellList->GetNumberOfIds();
|
|
|
|
if (nSendCells == 0)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
for (cellId = 0; cellId < nSendCells; cellId++)
|
|
{
|
|
subGridCellIds.insert(cellList->GetId(cellId));
|
|
}
|
|
}
|
|
|
|
globalPtIds[proc]->Delete();
|
|
|
|
int numUniqueCellIds = subGridCellIds.size();
|
|
|
|
if (numUniqueCellIds == 0)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
sendCells[proc]->SetNumberOfIds(numUniqueCellIds);
|
|
vtkIdType next = 0;
|
|
|
|
vtkstd::set<vtkIdType>::iterator it;
|
|
|
|
for (it = subGridCellIds.begin(); it != subGridCellIds.end(); ++it)
|
|
{
|
|
sendCells[proc]->SetId(next++, *it);
|
|
}
|
|
}
|
|
|
|
delete [] globalPtIds;
|
|
|
|
cellList->Delete();
|
|
|
|
return sendCells;
|
|
}
|
|
void vtkDistributedDataFilter::RemoveRemoteCellsFromList(
|
|
vtkIdList *cellList, int *gidCells, int *remoteCells, int nRemoteCells)
|
|
{
|
|
vtkIdType id, nextId;
|
|
int id2;
|
|
vtkIdType nLocalCells = cellList->GetNumberOfIds();
|
|
|
|
// both lists should be very small, so we just do an n^2 lookup
|
|
|
|
for (id = 0, nextId = 0; id < nLocalCells; id++)
|
|
{
|
|
vtkIdType localCellId = cellList->GetId(id);
|
|
int globalCellId = gidCells[localCellId];
|
|
|
|
int found = 0;
|
|
|
|
for (id2 = 0; id2 < nRemoteCells; id2++)
|
|
{
|
|
if (remoteCells[id2] == globalCellId)
|
|
{
|
|
found = 1;
|
|
break;
|
|
}
|
|
}
|
|
|
|
if (!found)
|
|
{
|
|
cellList->SetId(nextId++, localCellId);
|
|
}
|
|
}
|
|
|
|
cellList->SetNumberOfIds(nextId);
|
|
}
|
|
|
|
// Set the ghost levels for the points and cells in the received cells.
|
|
// Merge the new ghost cells into the supplied grid, and return the new grid.
|
|
// Delete all grids except the new merged grid.
|
|
|
|
vtkUnstructuredGrid *vtkDistributedDataFilter::SetMergeGhostGrid(
|
|
vtkUnstructuredGrid *ghostCellGrid,
|
|
vtkUnstructuredGrid *incomingGhostCells,
|
|
int ghostLevel, vtkDistributedDataFilterSTLCloak *idMap)
|
|
|
|
{
|
|
int i;
|
|
|
|
if (incomingGhostCells->GetNumberOfCells() < 1)
|
|
{
|
|
return ghostCellGrid;
|
|
}
|
|
|
|
// Set the ghost level of all new cells, and set the ghost level of all
|
|
// the points. We know some points in the new grids actually have ghost
|
|
// level one lower, because they were on the boundary of the previous
|
|
// grid. This is OK if ghostLevel is > 1. When we merge, vtkMergeCells
|
|
// will skip these points because they are already in the previous grid.
|
|
// But if ghostLevel is 1, those boundary points were in our original
|
|
// grid, and we need to use the global ID map to determine if the
|
|
// point ghost levels should be set to 0.
|
|
|
|
vtkDataArray *da = incomingGhostCells->GetCellData()->GetArray("vtkGhostLevels");
|
|
vtkUnsignedCharArray *cellGL = vtkUnsignedCharArray::SafeDownCast(da);
|
|
|
|
da = incomingGhostCells->GetPointData()->GetArray("vtkGhostLevels");
|
|
vtkUnsignedCharArray *ptGL = vtkUnsignedCharArray::SafeDownCast(da);
|
|
|
|
unsigned char *ia = cellGL->GetPointer(0);
|
|
|
|
for (i=0; i < incomingGhostCells->GetNumberOfCells(); i++)
|
|
{
|
|
ia[i] = (unsigned char)ghostLevel;
|
|
}
|
|
|
|
ia = ptGL->GetPointer(0);
|
|
|
|
for (i=0; i < incomingGhostCells->GetNumberOfPoints(); i++)
|
|
{
|
|
ia[i] = (unsigned char)ghostLevel;
|
|
}
|
|
|
|
// now merge
|
|
|
|
vtkUnstructuredGrid *mergedGrid = incomingGhostCells;
|
|
|
|
if (ghostCellGrid && (ghostCellGrid->GetNumberOfCells() > 0))
|
|
{
|
|
vtkDataSet *sets[2];
|
|
|
|
sets[0] = ghostCellGrid; // both sets will be deleted by MergeGrids
|
|
sets[1] = incomingGhostCells;
|
|
|
|
const char *nodeIds = this->GetGlobalNodeIdArrayName(ghostCellGrid);
|
|
|
|
mergedGrid =
|
|
vtkDistributedDataFilter::MergeGrids(sets, 2, DeleteYes, nodeIds, 0.0, NULL);
|
|
}
|
|
|
|
// If this is ghost level 1, mark any points from our original grid
|
|
// as ghost level 0.
|
|
|
|
if (ghostLevel == 1)
|
|
{
|
|
da = mergedGrid->GetPointData()->GetArray("vtkGhostLevels");
|
|
ptGL = vtkUnsignedCharArray::SafeDownCast(da);
|
|
|
|
int *gidPoints = this->GetGlobalNodeIds(mergedGrid);
|
|
int npoints = mergedGrid->GetNumberOfPoints();
|
|
|
|
vtkstd::map<int, int>::iterator imap;
|
|
|
|
for (i=0; i < npoints; i++)
|
|
{
|
|
imap = idMap->IntMap.find(gidPoints[i]);
|
|
|
|
if (imap != idMap->IntMap.end())
|
|
{
|
|
ptGL->SetValue(i,0); // found among my ghost level 0 cells
|
|
}
|
|
}
|
|
}
|
|
|
|
return mergedGrid;
|
|
}
|
|
vtkUnstructuredGrid *vtkDistributedDataFilter::MergeGrids(
|
|
vtkDataSet **sets, int nsets, int deleteDataSets,
|
|
const char *globalNodeIdArrayName, float pointMergeTolerance,
|
|
const char *globalCellIdArrayName)
|
|
{
|
|
int i;
|
|
|
|
if (nsets == 0)
|
|
{
|
|
return NULL;
|
|
}
|
|
|
|
vtkModelMetadata *mmd = NULL;
|
|
|
|
for (i=0; i<nsets; i++)
|
|
{
|
|
// It's possible we're merging regular cells (with metadata) with
|
|
// ghost cells (no metadata) so check each dataset for metadata.
|
|
|
|
int inputHasMetadata = vtkDistributedDataFilter::HasMetadata(sets[i]);
|
|
|
|
if (!inputHasMetadata)
|
|
{
|
|
continue;
|
|
}
|
|
|
|
vtkModelMetadata *submmd = vtkModelMetadata::New();
|
|
submmd->Unpack(sets[i], DeleteYes);
|
|
|
|
if (mmd)
|
|
{
|
|
mmd->MergeModelMetadata(submmd);
|
|
submmd->Delete();
|
|
}
|
|
else
|
|
{
|
|
mmd = submmd;
|
|
}
|
|
}
|
|
|
|
vtkUnstructuredGrid *newGrid = vtkUnstructuredGrid::New();
|
|
|
|
vtkMergeCells *mc = vtkMergeCells::New();
|
|
mc->SetUnstructuredGrid(newGrid);
|
|
|
|
mc->SetTotalNumberOfDataSets(nsets);
|
|
|
|
int totalPoints = 0;
|
|
int totalCells = 0;
|
|
|
|
for (i=0; i<nsets; i++)
|
|
{
|
|
totalPoints += sets[i]->GetNumberOfPoints();
|
|
totalCells += sets[i]->GetNumberOfCells();
|
|
}
|
|
|
|
mc->SetTotalNumberOfPoints(totalPoints);
|
|
mc->SetTotalNumberOfCells(totalCells);
|
|
|
|
if (globalNodeIdArrayName)
|
|
{
|
|
mc->SetGlobalIdArrayName(globalNodeIdArrayName);
|
|
}
|
|
else
|
|
{
|
|
mc->SetPointMergeTolerance(pointMergeTolerance);
|
|
}
|
|
|
|
if (globalCellIdArrayName)
|
|
{
|
|
mc->SetGlobalCellIdArrayName(globalCellIdArrayName);
|
|
}
|
|
|
|
for (i=0; i<nsets; i++)
|
|
{
|
|
mc->MergeDataSet(sets[i]);
|
|
|
|
if (deleteDataSets)
|
|
{
|
|
sets[i]->Delete();
|
|
}
|
|
}
|
|
|
|
mc->Finish();
|
|
mc->Delete();
|
|
|
|
if (mmd)
|
|
{
|
|
// Pack the metadata onto the new grid and delete it.
|
|
|
|
mmd->Pack(newGrid);
|
|
mmd->Delete();
|
|
}
|
|
|
|
return newGrid;
|
|
}
|
|
int vtkDistributedDataFilter::HasMetadata(vtkDataSet *s)
|
|
{
|
|
return vtkModelMetadata::HasMetadata(vtkUnstructuredGrid::SafeDownCast(s));
|
|
}
|
|
//-------------------------------------------------------------------------
|
|
int vtkDistributedDataFilter::FillInputPortInformation(int, vtkInformation *info)
|
|
{
|
|
info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkDataSet");
|
|
return 1;
|
|
}
|
|
|
|
//-------------------------------------------------------------------------
|
|
void vtkDistributedDataFilter::PrintSelf(ostream& os, vtkIndent indent)
|
|
{
|
|
this->Superclass::PrintSelf(os,indent);
|
|
|
|
os << indent << "Kdtree: " << this->Kdtree << endl;
|
|
os << indent << "Controller: " << this->Controller << endl;
|
|
os << indent << "NumProcesses: " << this->NumProcesses << endl;
|
|
os << indent << "MyId: " << this->MyId << endl;
|
|
os << indent << "Target: " << this->Target << endl;
|
|
os << indent << "Source: " << this->Source << endl;
|
|
if (this->GlobalNodeIdArrayName)
|
|
{
|
|
os << indent << "GlobalNodeIdArrayName: " << this->GlobalNodeIdArrayName << endl;
|
|
}
|
|
if (this->GlobalElementIdArrayName)
|
|
{
|
|
os << indent << "GlobalElementIdArrayName: " << this->GlobalElementIdArrayName << endl;
|
|
}
|
|
os << indent << "RetainKdtree: " << this->RetainKdtree << endl;
|
|
os << indent << "IncludeAllIntersectingCells: " << this->IncludeAllIntersectingCells << endl;
|
|
os << indent << "ClipCells: " << this->ClipCells << endl;
|
|
|
|
os << indent << "Timing: " << this->Timing << endl;
|
|
os << indent << "UseMinimalMemory: " << this->UseMinimalMemory << endl;
|
|
}
|
|
|
|
|