src/Plugins/GhostCellsGenerator/vtkPUnstructuredGridGhostCellsGenerator.cxx

   1 /*=========================================================================
   2
   3   Program:   Visualization Toolkit
   4   Module:    vtkPUnstructuredGridGhostCellsGenerator.cxx
   5
   6   Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
   7   All rights reserved.
   8   See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
   9
  10      This software is distributed WITHOUT ANY WARRANTY; without even
  11      the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  12      PURPOSE.  See the above copyright notice for more information.
  13
  14 =========================================================================*/
  15
  16 #include "vtkPUnstructuredGridGhostCellsGenerator.h"
  17
  18 #include "vtkCellArray.h"
  19 #include "vtkCharArray.h"
  20 #include "vtkDataSetSurfaceFilter.h"
  21 #include "vtkExtractCells.h"
  22 #include "vtkIdList.h"
  23 #include "vtkIdTypeArray.h"
  24 #include "vtkInformation.h"
  25 #include "vtkInformationVector.h"
  26 #include "vtkMergeCells.h"
  27 #include "vtkMergePoints.h"
  28 #include "vtkMPICommunicator.h"
  29 #include "vtkMultiProcessController.h"
  30 #include "vtkNew.h"
  31 #include "vtkObjectFactory.h"
  32 #include "vtkPointData.h"
  33 #include "vtkPoints.h"
  34 #include "vtkStreamingDemandDrivenPipeline.h"
  35 #include "vtkUnstructuredGrid.h"
  36
  37 #include <map>
  38 #include <set>
  39 #include <vector>
  40
  41 //----------------------------------------------------------------------------
  42 // Helpers
  43 namespace
  44 {
  45 template<class T>
  46 bool AllGatherV(vtkMultiProcessController* controller,
  47                 const T* localV,
  48                 vtkIdType localSize,
  49                 std::vector<T>& globalV,
  50                 std::vector<vtkIdType>& sizes,
  51                 std::vector<vtkIdType>& offsets)
  52 {
  53   int nbOfRanks = controller->GetNumberOfProcesses();
  54   sizes.resize(nbOfRanks);
  55   int ret = controller->AllGather(&localSize, &sizes[0], 1);
  56   if (ret == 0)
  57     {
  58     vtkErrorWithObjectMacro(controller, << "Communication error!");
  59     return false;
  60     }
  61   vtkIdType count = 0;
  62   offsets.resize(nbOfRanks);
  63   for (int i = 0; i < nbOfRanks; i++)
  64     {
  65     offsets[i] = count;
  66     count += sizes[i];
  67     }
  68   globalV.resize(count);
  69   if (count > 0)
  70     {
  71     controller->AllGatherV(localSize > 0 ? localV : 0,
  72       &globalV[0], localSize, &sizes[0], &offsets[0]);
  73     }
  74   return true;
  75 }
  76 }
  77
  78 //----------------------------------------------------------------------------
  79 // Internal data structures
  80
  81 // Class to hold asynchronous communication information
  82 class CommDataInfo
  83 {
  84 public:
  85   CommDataInfo() : SendLen(-1), RecvLen(-1), CommStep(0)
  86   {
  87     this->SendBuffer = vtkCharArray::New();
  88     this->RecvBuffer = vtkCharArray::New();
  89   }
  90
  91   CommDataInfo(const CommDataInfo& c)
  92   {
  93     *this = c;
  94     if (this->SendBuffer) { this->SendBuffer->Register(0); }
  95     if (this->RecvBuffer) { this->RecvBuffer->Register(0); }
  96   }
  97
  98   ~CommDataInfo()
  99   {
 100     if (this->SendBuffer) { this->SendBuffer->Delete(); }
 101     if (this->RecvBuffer) { this->RecvBuffer->Delete(); }
 102   }
 103
 104   vtkMPICommunicator::Request SendReqs[2];
 105   vtkMPICommunicator::Request RecvReqs[2];
 106   vtkCharArray *SendBuffer;
 107   vtkCharArray *RecvBuffer;
 108   vtkIdType SendLen;
 109   vtkIdType RecvLen;
 110   int CommStep;
 111 };
 112
 113 // Communication arrays
 114 struct vtkPUnstructuredGridGhostCellsGenerator::vtkInternals
 115 {
 116   // For global ids
 117   std::map<vtkIdType, vtkIdType> GlobalToLocalPointIdMap;
 118   std::vector<vtkIdType> AllGlobalIdsOfSurfacePoints;
 119
 120   // For point coordinates
 121   vtkNew<vtkMergePoints> LocalPoints;
 122   std::vector<vtkIdType> LocalPointsMap;
 123   std::vector<double> AllPointsOfSurfacePoints;
 124
 125   std::vector<vtkIdType> AllSizes;
 126   std::vector<vtkIdType> AllOffsets;
 127   std::map<int, std::set<vtkIdType> > NeighborRanksCells;
 128   std::map<int, CommDataInfo> CommData;
 129   vtkUnstructuredGridBase* Input;
 130
 131   vtkDataArray* InputGlobalPointIds;
 132   bool UseGlobalPointIds;
 133 };
 134
 135 static const int UGGCG_SIZE_EXCHANGE_TAG = 9000;
 136 static const int UGGCG_DATA_EXCHANGE_TAG = 9001;
 137 static const char* UGGCG_GLOBAL_POINT_IDS = "GlobalNodeIds";
 138
 139 //----------------------------------------------------------------------------
 140
 141 vtkStandardNewMacro(vtkPUnstructuredGridGhostCellsGenerator)
 142 vtkSetObjectImplementationMacro(
 143   vtkPUnstructuredGridGhostCellsGenerator, Controller, vtkMultiProcessController);
 144
 145 //----------------------------------------------------------------------------
 146 vtkPUnstructuredGridGhostCellsGenerator::vtkPUnstructuredGridGhostCellsGenerator()
 147 {
 148   this->Controller = NULL;
 149   this->SetController(vtkMultiProcessController::GetGlobalController());
 150
 151   this->Internals = NULL;
 152   this->GlobalPointIdsArrayName = NULL;
 153   this->UseGlobalPointIds = true;
 154   this->SetGlobalPointIdsArrayName(UGGCG_GLOBAL_POINT_IDS);
 155 }
 156
 157 //----------------------------------------------------------------------------
 158 vtkPUnstructuredGridGhostCellsGenerator::~vtkPUnstructuredGridGhostCellsGenerator()
 159 {
 160   this->SetController(NULL);
 161   this->SetGlobalPointIdsArrayName(NULL);
 162
 163   delete this->Internals;
 164   this->Internals = 0;
 165 }
 166
 167 //-----------------------------------------------------------------------------
 168 void vtkPUnstructuredGridGhostCellsGenerator::PrintSelf(ostream& os, vtkIndent indent)
 169 {
 170   Superclass::PrintSelf(os, indent);
 171 }
 172
 173 //-----------------------------------------------------------------------------
 174 int vtkPUnstructuredGridGhostCellsGenerator::RequestData(
 175   vtkInformation *vtkNotUsed(request),
 176   vtkInformationVector **inputVector,
 177   vtkInformationVector *outputVector)
 178 {
 179   // get the info objects
 180   vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
 181   vtkInformation *outInfo = outputVector->GetInformationObject(0);
 182
 183   // get the input and output. Input may just have the UnstructuredGridBase
 184   // interface, but output should be an unstructured grid.
 185   vtkUnstructuredGridBase *input = vtkUnstructuredGridBase::SafeDownCast(
 186     inInfo->Get(vtkDataObject::DATA_OBJECT()));
 187   vtkUnstructuredGrid *output = vtkUnstructuredGrid::SafeDownCast(
 188     outInfo->Get(vtkDataObject::DATA_OBJECT()));
 189
 190   int ghostLevels = 1;
 191   //ghostLevels = outInfo->Get(
 192   //  vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS());
 193
 194   if (!this->Controller)
 195     {
 196     this->Controller = vtkMultiProcessController::GetGlobalController();
 197     }
 198   output->Reset();
 199
 200   this->NumRanks = this->Controller ? this->Controller->GetNumberOfProcesses() : 1;
 201   this->RankId = this->Controller ? this->Controller->GetLocalProcessId() : 0;
 202   if (ghostLevels == 0 || !this->Controller || this->NumRanks == 1)
 203     {
 204     // Ghost levels are not requested. Nothing to do but pass the dataset.
 205     if (this->RankId == 0)
 206       {
 207       vtkWarningMacro(<< "Ghost cells not requested or not needed.");
 208       }
 209     output->ShallowCopy(input);
 210     return 1;
 211     }
 212
 213   delete this->Internals;
 214   this->Internals = new vtkPUnstructuredGridGhostCellsGenerator::vtkInternals();
 215   this->Internals->Input = input;
 216
 217   vtkPointData *inputPD = input->GetPointData();
 218   this->Internals->InputGlobalPointIds = inputPD->GetGlobalIds();
 219
 220   if (!this->Internals->InputGlobalPointIds)
 221     {
 222     this->Internals->InputGlobalPointIds =
 223       inputPD->GetArray(this->GlobalPointIdsArrayName);
 224     inputPD->SetGlobalIds(this->Internals->InputGlobalPointIds);
 225     }
 226
 227   if (!this->UseGlobalPointIds)
 228     {
 229     this->Internals->InputGlobalPointIds = NULL;
 230     }
 231
 232   this->ExtractAndReduceSurfacePoints();
 233   this->UpdateProgress(0.3);
 234
 235   this->ComputeSharedPoints();
 236   this->UpdateProgress(0.6);
 237
 238   this->ExtractAndSendGhostCells();
 239   this->UpdateProgress(0.8);
 240
 241   this->ReceiveAndMergeGhostCells(output);
 242   this->UpdateProgress(1.0);
 243
 244   output->GetInformation()->Set(vtkDataObject::DATA_NUMBER_OF_GHOST_LEVELS(), 1);
 245
 246   this->Controller->Barrier();
 247
 248   delete this->Internals;
 249   this->Internals = 0;
 250
 251   return 1;
 252 }
 253
 254 //-----------------------------------------------------------------------------
 255 // Step 1: Extract surface geometry and all reduce global ids of surface points
 256 void vtkPUnstructuredGridGhostCellsGenerator::ExtractAndReduceSurfacePoints()
 257 {
 258   // Extract boundary cells and points with the surface filter
 259   vtkNew<vtkDataSetSurfaceFilter> surfaceFilter;
 260   surfaceFilter->SetInputData(this->Internals->Input);
 261   surfaceFilter->PassThroughPointIdsOn();
 262   surfaceFilter->Update();
 263
 264   vtkPolyData* surface = surfaceFilter->GetOutput();
 265   vtkIdType nbSurfacePoints = surface->GetNumberOfPoints();
 266   vtkCellArray* surfaceCells = surface->GetPolys();
 267   surfaceCells->InitTraversal();
 268   vtkIdType npts, *pts;
 269
 270   vtkIdTypeArray* surfaceOriginalPointIds = vtkIdTypeArray::SafeDownCast(
 271     surface->GetPointData()->GetArray(surfaceFilter->GetOriginalPointIdsName()));
 272
 273   if (this->Internals->InputGlobalPointIds)
 274     {
 275     std::vector<vtkIdType> globalIdsOfSurfacePoints;
 276     globalIdsOfSurfacePoints.reserve(nbSurfacePoints);
 277
 278     // Browse surface cells and save global and local ids of cell points
 279     while (surfaceCells->GetNextCell(npts, pts))
 280       {
 281       for (vtkIdType i = 0; i < npts; i++)
 282         {
 283         vtkIdType origPtId = surfaceOriginalPointIds->GetValue(pts[i]);
 284         vtkIdType globalPtId = static_cast<vtkIdType>(
 285           this->Internals->InputGlobalPointIds->GetTuple1(origPtId));
 286
 287         if (this->Internals->GlobalToLocalPointIdMap.find(globalPtId) ==
 288           this->Internals->GlobalToLocalPointIdMap.end())
 289           {
 290           this->Internals->GlobalToLocalPointIdMap[globalPtId] = origPtId;
 291           globalIdsOfSurfacePoints.push_back(globalPtId);
 292           }
 293         }
 294       }
 295
 296     // Now reduce surface point global ids on ALL ranks
 297     ::AllGatherV(this->Controller, &globalIdsOfSurfacePoints[0],
 298       globalIdsOfSurfacePoints.size(),
 299       this->Internals->AllGlobalIdsOfSurfacePoints,
 300       this->Internals->AllSizes, this->Internals->AllOffsets);
 301     }
 302   else
 303     {
 304     // We can't use global ids, so we will process point coordinates instead
 305     vtkPoints* inputPoints = this->Internals->Input->GetPoints();
 306     vtkNew<vtkPoints> surfacePoints;
 307     surfacePoints->SetDataTypeToDouble();
 308     surfacePoints->Allocate(nbSurfacePoints);
 309     this->Internals->LocalPoints->InitPointInsertion(
 310       surfacePoints.Get(), surface->GetPoints()->GetBounds());
 311     this->Internals->LocalPointsMap.reserve(nbSurfacePoints);
 312
 313     // Browse surface cells and push point coordinates to the locator
 314     while (surfaceCells->GetNextCell(npts, pts))
 315       {
 316       for (vtkIdType i = 0; i < npts; i++)
 317         {
 318         vtkIdType origPtId = surfaceOriginalPointIds->GetValue(pts[i]);
 319         double p[3];
 320         inputPoints->GetPoint(origPtId, p);
 321         vtkIdType sid;
 322         if (this->Internals->LocalPoints->InsertUniquePoint(p, sid))
 323           {
 324           // New point, save the id of the original grid point id associated
 325           // to this surface point
 326           if (static_cast<vtkIdType>(this->Internals->LocalPointsMap.size()) <= sid)
 327             {
 328             this->Internals->LocalPointsMap.resize(sid + 1);
 329             }
 330           this->Internals->LocalPointsMap[sid] = origPtId;
 331           }
 332         }
 333       }
 334
 335     // Now reduce surface point coordinates on ALL ranks
 336     ::AllGatherV(this->Controller,
 337       (double*)surfacePoints->GetVoidPointer(0),
 338       surfacePoints->GetNumberOfPoints() * 3,
 339       this->Internals->AllPointsOfSurfacePoints,
 340       this->Internals->AllSizes, this->Internals->AllOffsets);
 341     }
 342 }
 343
 344 //---------------------------------------------------------------------------
 345 // Step 2: browse global ids/point coordinates of other ranks and check if some
 346 // are duplicated locally.
 347 // For each neighbor rank, save the ids of the cells adjacent to the surface
 348 // points shared, those cells are the ghost cells we will send them.
 349 void vtkPUnstructuredGridGhostCellsGenerator::ComputeSharedPoints()
 350 {
 351   vtkNew<vtkIdList> cellIdsList;
 352   for (int i = 0; i < this->NumRanks; i++)
 353     {
 354     if (i == this->RankId)
 355       {
 356       continue;
 357       }
 358     for (vtkIdType j = 0, idx = this->Internals->AllOffsets[i];
 359          j < this->Internals->AllSizes[i]; j++, idx++)
 360       {
 361       vtkIdType localPointId = -1;
 362       if (this->Internals->InputGlobalPointIds)
 363         {
 364         // Check if this point exists locally from its global ids, if so
 365         // get its local id.
 366         vtkIdType gid = this->Internals->AllGlobalIdsOfSurfacePoints[idx];
 367         std::map<vtkIdType, vtkIdType>::iterator iter =
 368           this->Internals->GlobalToLocalPointIdMap.find(gid);
 369         if (iter != this->Internals->GlobalToLocalPointIdMap.end())
 370           {
 371           localPointId = iter->second;
 372           }
 373         }
 374       else
 375         {
 376         // Check if this point exists locally from its coordinates, if so
 377         // get its local id.
 378         double *p = &this->Internals->AllPointsOfSurfacePoints[idx];
 379         localPointId = this->Internals->LocalPoints->IsInsertedPoint(p);
 380
 381         if (localPointId != -1)
 382           {
 383           localPointId = this->Internals->LocalPointsMap[localPointId];
 384           }
 385         idx += 2; // jump to next coordinates
 386         j += 2;
 387         }
 388
 389       if (localPointId != -1)
 390         {
 391         // Current rank also has a copy of this global point
 392         cellIdsList->Reset();
 393         // Get the cells connected to this point
 394         this->Internals->Input->GetPointCells(localPointId, cellIdsList.Get());
 395         vtkIdType nbIds = cellIdsList->GetNumberOfIds();
 396         // Add those cells to the list of cells to send to this rank
 397         for (vtkIdType k = 0; k < nbIds; k++)
 398           {
 399           this->Internals->NeighborRanksCells[i].insert(cellIdsList->GetId(k));
 400           }
 401         }
 402       }
 403     }
 404
 405   // Release memory of all reduced arrays
 406   this->Internals->AllGlobalIdsOfSurfacePoints.resize(0);
 407   this->Internals->AllPointsOfSurfacePoints.resize(0);
 408   this->Internals->AllSizes.resize(0);
 409   this->Internals->AllOffsets.resize(0);
 410   // Now we know our neighbors and which points we have in common and the
 411   // ghost cells to share.
 412 }
 413
 414 //-----------------------------------------------------------------------------
 415 // Step 3: extract and send the ghost cells to the neighbor ranks
 416 void vtkPUnstructuredGridGhostCellsGenerator::ExtractAndSendGhostCells()
 417 {
 418   vtkNew<vtkIdList> cellIdsList;
 419   vtkNew<vtkExtractCells> extractCells;
 420   extractCells->SetInputData(this->Internals->Input);
 421   std::map<int, std::set<vtkIdType> >::iterator iter =
 422     this->Internals->NeighborRanksCells.begin();
 423
 424   vtkMPICommunicator* com =
 425     vtkMPICommunicator::SafeDownCast(this->Controller->GetCommunicator());
 426
 427   // Browse all neighbor ranks and extract the cells connected to the points we share
 428   for (; iter != this->Internals->NeighborRanksCells.end(); ++iter)
 429     {
 430     int toRank = iter->first;
 431     std::set<vtkIdType>& cellsToShare = iter->second;
 432     cellIdsList->SetNumberOfIds(cellsToShare.size());
 433     std::set<vtkIdType>::iterator sIter = cellsToShare.begin();
 434     for (vtkIdType i = 0; sIter != cellsToShare.end(); ++sIter, i++)
 435       {
 436       cellIdsList->SetId(i, *sIter);
 437       }
 438     extractCells->SetCellList(cellIdsList.Get());
 439     extractCells->Update();
 440     vtkUnstructuredGrid* extractGrid = extractCells->GetOutput();
 441
 442     // Send the extracted grid to the neighbor rank asynchronously
 443     CommDataInfo& c = this->Internals->CommData[toRank];
 444     if (vtkCommunicator::MarshalDataObject(extractGrid, c.SendBuffer))
 445       {
 446       c.SendLen = c.SendBuffer->GetNumberOfTuples();
 447       // Send data length
 448       com->NoBlockSend(&c.SendLen, 1, toRank, UGGCG_SIZE_EXCHANGE_TAG, c.SendReqs[0]);
 449       // Send raw data
 450       com->NoBlockSend((char*)c.SendBuffer->GetVoidPointer(0), c.SendLen, toRank,
 451         UGGCG_DATA_EXCHANGE_TAG, c.SendReqs[1]);
 452       }
 453     }
 454 }
 455
 456 //-----------------------------------------------------------------------------
 457 // Step 4: receive the ghost cells from the neighbor ranks and merge them
 458 // to the local grid.
 459 void vtkPUnstructuredGridGhostCellsGenerator::ReceiveAndMergeGhostCells(
 460   vtkUnstructuredGrid *output)
 461 {
 462   // We need to compute a rough estimation of the total number of cells and
 463   // points for vtkMergeCells
 464   vtkIdType totalNbCells = this->Internals->Input->GetNumberOfCells();
 465   vtkIdType totalNbPoints = this->Internals->Input->GetNumberOfPoints();
 466
 467   vtkMPICommunicator* com =
 468     vtkMPICommunicator::SafeDownCast(this->Controller->GetCommunicator());
 469
 470   // Browse all neighbor ranks and receive the mesh that contains cells
 471   int nbNeighbors = static_cast<int>(this->Internals->NeighborRanksCells.size());
 472   std::vector<vtkUnstructuredGridBase*> neighborGrids;
 473   neighborGrids.reserve(nbNeighbors);
 474
 475   // First create requests to receive the size of the mesh to receive
 476   std::map<int, std::set<vtkIdType> >::iterator iter;
 477   for (iter = this->Internals->NeighborRanksCells.begin();
 478     iter != this->Internals->NeighborRanksCells.end(); ++iter)
 479     {
 480     vtkIdType fromRank = iter->first;
 481     CommDataInfo& c = this->Internals->CommData[fromRank];
 482     com->NoBlockReceive(
 483       &c.RecvLen, 1, fromRank, UGGCG_SIZE_EXCHANGE_TAG, c.RecvReqs[0]);
 484     }
 485
 486   // Then, once the data length is received, create requests to receive the mesh data
 487   int counter = 0;
 488   while (counter != nbNeighbors)
 489     {
 490     for (iter = this->Internals->NeighborRanksCells.begin();
 491       iter != this->Internals->NeighborRanksCells.end(); ++iter)
 492       {
 493       vtkIdType fromRank = iter->first;
 494       CommDataInfo& c = this->Internals->CommData[fromRank];
 495       if (!c.RecvReqs[0].Test() || c.CommStep != 0)
 496         {
 497         continue;
 498         }
 499       c.CommStep = 1;
 500       c.RecvBuffer->SetNumberOfValues(c.RecvLen);
 501       com->NoBlockReceive(
 502         (char*)c.RecvBuffer->GetVoidPointer(0), c.RecvLen, fromRank,
 503         UGGCG_DATA_EXCHANGE_TAG, c.RecvReqs[1]);
 504       counter++;
 505       }
 506     }
 507
 508   // Browse all neighbor ranks and receive the mesh that contains cells
 509   // that are ghost cells for current rank.
 510   counter = 0;
 511   while (counter != nbNeighbors)
 512     {
 513     for (iter = this->Internals->NeighborRanksCells.begin();
 514     iter != this->Internals->NeighborRanksCells.end(); ++iter)
 515       {
 516       vtkIdType fromRank = iter->first;
 517       CommDataInfo& c = this->Internals->CommData[fromRank];
 518
 519       if (!c.RecvReqs[1].Test() || c.CommStep != 1)
 520         {
 521         continue;
 522         }
 523
 524       c.CommStep = 2;
 525       vtkUnstructuredGrid* grid = vtkUnstructuredGrid::New();
 526       vtkCommunicator::UnMarshalDataObject(c.RecvBuffer, grid);
 527       c.RecvBuffer->Delete();
 528       c.RecvBuffer = NULL;
 529
 530       // Flag the received grid elements as ghosts
 531       grid->AllocatePointGhostArray();
 532       grid->AllocateCellGhostArray();
 533       grid->GetPointGhostArray()->FillComponent(0, 1);
 534       grid->GetCellGhostArray()->FillComponent(0, 1);
 535
 536       // Make sure the global point ids array is tagged accordingly
 537       if (this->Internals->InputGlobalPointIds &&
 538           !grid->GetPointData()->GetGlobalIds())
 539         {
 540         grid->GetPointData()->SetGlobalIds(grid->GetPointData()->GetArray(
 541           this->Internals->InputGlobalPointIds->GetName()));
 542         }
 543
 544       totalNbCells += grid->GetNumberOfCells();
 545       totalNbPoints += grid->GetNumberOfPoints();
 546
 547       neighborGrids.push_back(grid);
 548
 549       counter++;
 550       }
 551     }
 552
 553   // Shallow copy the input grid and initialize the ghost arrays
 554   vtkNew<vtkUnstructuredGrid> inputCopy;
 555   inputCopy->ShallowCopy(this->Internals->Input);
 556   inputCopy->AllocatePointGhostArray();
 557   inputCopy->AllocateCellGhostArray();
 558
 559   // MergeCells merge input + grids that contains ghost cells to the output grid
 560   vtkNew<vtkMergeCells> mergeCells;
 561   mergeCells->SetUnstructuredGrid(output);
 562   mergeCells->SetTotalNumberOfCells(totalNbCells);
 563   mergeCells->SetTotalNumberOfPoints(totalNbPoints);
 564   mergeCells->SetTotalNumberOfDataSets(
 565     1 + static_cast<int>(this->Internals->NeighborRanksCells.size()));
 566   mergeCells->SetUseGlobalIds(this->Internals->InputGlobalPointIds != 0 ? 1 : 0);
 567   mergeCells->SetPointMergeTolerance(0.0);
 568
 569   // Merge input grid first
 570   mergeCells->MergeDataSet(inputCopy.Get());
 571
 572   // Then merge ghost grid from neighbor rank
 573   for (std::size_t i = 0; i < neighborGrids.size(); i++)
 574     {
 575     mergeCells->MergeDataSet(neighborGrids[i]);
 576     neighborGrids[i]->Delete();
 577     }
 578
 579   // Finalize the merged output
 580   mergeCells->Finish();
 581 }