feat: Add new scatter strategies for initial mesh redistribution by bd713 · Pull Request #4033 · GEOS-DEV/GEOS

bd713 · 2026-04-21T03:34:53Z

Problem

The current strategy to scatter the mesh initially loaded by rank 0 uses vtkRedistributeDataSetFilter (kd-tree partitioning). This approach can take several minutes for large meshes (10s of millions of elements) when redistributed among 1000+ ranks.

Recent example: A 136M cell mesh redistributed to 1024 ranks takes about 17 minutes with the current kd-tree approach.

Solution

This PR adds partitioning strategies that bypass the VTK library's bottleneck, available via the new scatterMethod XML attribute.

Available strategies

Method	Description	Comment
`kdtree`	VTK built-in kd-tree (default)	Spatial locality, but can be slow and current version can lose tetra or pyramid cells
`contiguous`	Contiguous cell-index blocks	Fast, but does not guarantee spatial locality
`cartesian`	Regular Cartesian grid using `-x/-y/-z` partition counts	Does not guarantee load balancing if mesh is irregular
`rcb`	Recursive Coordinate Bisection along longest bounding-box axis	Spatial locality + load balance

Remark

The resulting partitions can be further refined and improved via calls to ParMetis or PTScotch with partitionRefinement>0

Benchmark

==============================================================================================================
SUMMARY (1024 ranks, 135901500 cells)
==============================================================================================================
Method        | Time (s)  | Cells       | Status  | Min cells | Max cells | Imbalance   | BBox Overlap
--------------------------------------------------------------------------------------------------------------
KdTree        | 1020.194  | 135901494   | LOSS    | 126714    | 139764    | 1.053103    | 1.086       
Contiguous    |   15.822  | 135901500   | OK      | 132716    | 132717    | 1.000005    | 2.861       
Cartesian     |   22.005  | 135901500   | OK      | 116157    | 152266    | 1.147304    | 1.095       
RCB           |   34.289  | 135901500   | OK      | 132716    | 132717    | 1.000005    | 1.116       
--------------------------------------------------------------------------------------------------------------

Speedup relative to each method:
  KdTree        64.48x slower than fastest
  Contiguous    1.00x slower than fastest <-- fastest
  Cartesian     1.39x slower than fastest
  RCB           2.17x slower than fastest

Best load balance: Contiguous (imbalance ratio: 1.000005x)
Best compactness:  KdTree (bbox overlap: 1.086)
==============================================================================================================

Rebaseline needed

67/78 failures are XML metadata diff: Group has a child 'scatterMethod' in the file to compare but not the baseline file.
11/78 are real partitioning changes in contact cases as we switched from Morton's algorithm to RCB

…tter

victorapm · 2026-04-23T05:32:15Z

That looks great, Bertrand. Even though contiguous is the fastest, we might want to have RCB as default since bbox values look much better and the additional cost is not too bad. cc @castelletto1

dkachuma

Looks good to me.

…tter

…GEOS into feat/bd713/meshScatter

joshua-white

Nice. Only suggestion might be to add a log note if a "large" mesh and kdtree are detected, suggesting user consider something else.

# Conflicts: # src/coreComponents/integrationTests/meshTests/testVTKImport.cpp # src/coreComponents/mesh/generators/VTKUtilities.cpp

…tter

victorapm · 2026-10-27T02:29:29Z

Needs code owner review from @OmarDuran , @jafranc and/or @rrsettgast Thanks 🙏

jafranc

Looks good, but some part need to be moved to honor responsability logic at header level a bit more (and avoid too much spreading across the code base)

jafranc · 2026-10-27T15:43:40Z

  // Identify which 3D global IDs we need from other ranks
-  SortedArray< int64_t > missingGlobalIds;
+  stdVector< int64_t > missingGlobalIdsVec;


Now that we are using c++20

auto missing = neighbors2Dto3D | std::views::join | std::views::filter( [&]( int64_t id ){ return !localPartitionMap.contains( id ); } ); stdVector< int64_t > missingGlobalIdsVec( missing.begin(), missing.end() ); // C++20 iterator-pair ctor std::ranges::sort( missingGlobalIdsVec ); auto const tail = std::ranges::unique( missingGlobalIdsVec ); missingGlobalIdsVec.erase( tail.begin(), tail.end() );

sounds simpler.

jafranc · 2026-10-27T16:54:42Z

+stdVector< stdArray< real64, 3 > >
+computeCentroids( vtkDataSet & mesh )
+{
+  vtkIdType const n = mesh.GetNumberOfCells();
+  stdVector< stdArray< real64, 3 > > centroids( n );
+  vtkNew< vtkIdList > ptIds;
+
+  for( vtkIdType c = 0; c < n; ++c )
+  {
+    mesh.GetCellPoints( c, ptIds );
+    real64 cx = 0.0, cy = 0.0, cz = 0.0;
+    vtkIdType const nPts = ptIds->GetNumberOfIds();
+    for( vtkIdType i = 0; i < nPts; ++i )
+    {
+      real64 p[3];
+      mesh.GetPoint( ptIds->GetId( i ), p );
+      cx += p[0];
+      cy += p[1];
+      cz += p[2];
+    }
+    if( nPts > 0 )
+    {
+      real64 const inv = 1.0 / nPts;
+      centroids[c] = { cx * inv, cy * inv, cz * inv };
+    }
+    else
+    {
+      centroids[c] = { 0.0, 0.0, 0.0 };
+    }
+  }
+  return centroids;
+}


can't we reuse something from ComputationalGeom here ? or as we are using VTK the canonical VTK function for that ?

As for now, centroids can fall outside the cells.

jafranc · 2026-10-28T07:56:53Z

+// Buffer serialization helpers
+// ============================================================================
+
+void appendBytes( stdVector< char > & buf, void const * data, int64_t n )


why won't these fit in mesh/BufferOps.hpp ?

jafranc · 2026-10-28T08:06:16Z

+// MPI large-message helpers (handles buffers > 2 GB)
+// ============================================================================
+
+constexpr int64_t MPI_CHUNK_SIZE = INT64_C( 1 ) << 30; // 1 GB


same here, why isn't it in MpiWrapper.hpp or even better directly a logic of MpiWrapper::send on the maximal chunck size and checkSum ?

jafranc · 2026-10-28T08:23:32Z

+  cells->ConvertToDefaultStorage();
+
+  vtkDataArray * offsets = cells->GetOffsetsArray();
+  vtkDataArray * conn = cells->GetConnectivityArray();


Idk the status elsewhere in the code and if we should guard this with VTK version warning and throw otherwise or not

jafranc · 2026-10-28T08:27:02Z

+    vtkSmartPointer< vtkCellArray > cellArray = readCellArray( ptr );
+    vtkSmartPointer< vtkCellArray > faceLocations = readCellArray( ptr );
+    vtkSmartPointer< vtkCellArray > faces = readCellArray( ptr );


Suggested change

vtkSmartPointer< vtkCellArray > cellArray = readCellArray( ptr );

vtkSmartPointer< vtkCellArray > faceLocations = readCellArray( ptr );

vtkSmartPointer< vtkCellArray > faces = readCellArray( ptr );

//sweep ptr reading

vtkSmartPointer< vtkCellArray > cellArray = readCellArray( ptr );

vtkSmartPointer< vtkCellArray > faceLocations = readCellArray( ptr );

vtkSmartPointer< vtkCellArray > faces = readCellArray( ptr );

I would add an extra reading sweep function as this is implicit convention and will be hard to maintain.

jafranc · 2026-10-28T10:46:02Z

+  for( vtkIdType i = 0; i < static_cast< vtkIdType >( assignment.size() ); ++i )
+  {
+    if( assignment[i] < mid )
+    {
+      loCells.push_back( i );
+      loAssign.push_back( assignment[i] );
+    }
+    else
+    {
+      hiCells.push_back( i );
+      hiAssign.push_back( assignment[i] );
+    }
+  }


Suggested change

for( vtkIdType i = 0; i < static_cast< vtkIdType >( assignment.size() ); ++i )

{

if( assignment[i] < mid )

{

loCells.push_back( i );

loAssign.push_back( assignment[i] );

}

else

{

hiCells.push_back( i );

hiAssign.push_back( assignment[i] );

}

}

auto const n = static_cast< vtkIdType >( assignment.size() );

stdVector< vtkIdType > order( n );

std::iota( order.begin(), order.end(), vtkIdType{ 0 } );

auto const split = std::stable_partition( order.begin(), order.end(),

[&]( vtkIdType i ) { return assignment[i] < mid; } );

loCells.assign( order.begin(), split );

hiCells.assign( split, order.end() );

auto const rankOf = [&]( vtkIdType i ) { return assignment[i]; };

std::ranges::transform( loCells, std::back_inserter( loAssign ), rankOf );

std::ranges::transform( hiCells, std::back_inserter( hiAssign ), rankOf );

jafranc · 2026-10-28T10:49:15Z

+  if( hiCells.size() <= loCells.size() )
+  {
+    result.hiMesh = extract( mesh, hiCells );
+    result.hiAssignment = std::move( hiAssign );
+    result.loMesh = extract( mesh, loCells );
+    result.loAssignment = std::move( loAssign );
+  }
+  else
+  {
+    result.loMesh = extract( mesh, loCells );
+    result.loAssignment = std::move( loAssign );
+    result.hiMesh = extract( mesh, hiCells );
+    result.hiAssignment = std::move( hiAssign );
+  }


Maybe move this to the struct ctor would be cleaner ?

jafranc · 2026-10-28T11:03:19Z

+  std::function< void( vtkIdType, vtkIdType, integer, integer ) > bisect;
+  bisect = [&]( vtkIdType begin, vtkIdType end, integer rankLo, integer rankHi )


Suggested change

std::function< void( vtkIdType, vtkIdType, integer, integer ) > bisect;

bisect = [&]( vtkIdType begin, vtkIdType end, integer rankLo, integer rankHi )

auto bisect = [&]( vtkIdType begin, vtkIdType end, integer rankLo, integer rankHi )

jafranc · 2026-10-28T11:24:44Z

+  // KdTree: legacy path using VTK's built-in redistribution
+  if( method == ScatterMethod::kdtree )
+  {
+    constexpr vtkIdType largeMeshCellCount = 10000000;


looks like it should be a capital const

DENEL Bertrand added 2 commits April 20, 2026 22:01

Draft

e7eae8c

Merge remote-tracking branch 'origin/develop' into feat/bd713/meshSca…

1549f0e

…tter

bd713 self-assigned this Apr 21, 2026

bd713 added the type: feature New feature or request label Apr 21, 2026

bd713 requested a review from victorapm April 21, 2026 03:35

bd713 added ci: run CUDA builds ci: run integrated tests Allows to run the integrated tests in GEOS CI labels Apr 27, 2026

bd713 marked this pull request as ready for review April 27, 2026 17:09

bd713 requested review from OmarDuran, jafranc, rrsettgast and wrtobin as code owners April 27, 2026 17:09

Fixed dependency

73a367e

bd713 added flag: requires rebaseline Requires rebaseline branch in integratedTests ci: run code coverage enables running of the code coverage CI jobs labels Apr 27, 2026

DENEL Bertrand added 2 commits May 20, 2026 16:18

merged

302f5c9

Merged with superCells

927b190

bd713 requested review from MelReyCG, corbett5 and dkachuma as code owners May 21, 2026 01:46

missing file

e7d3d73

dkachuma approved these changes May 22, 2026

View reviewed changes

bd713 and others added 4 commits July 22, 2026 20:03

Merge branch 'develop' into feat/bd713/meshScatter

379a3e0

Fix sorting and benchmark

534d027

Merge remote-tracking branch 'origin/develop' into feat/bd713/meshSca…

d5eef3a

…tter

Merge branch 'feat/bd713/meshScatter' of https://github.com/GEOS-DEV/…

7f98da1

…GEOS into feat/bd713/meshScatter

bd713 requested review from frankfeifan, herve-gross and jhuang2601 as code owners July 28, 2026 01:27

bd713 requested review from joshua-white and tjb-ltk as code owners July 28, 2026 01:27

feat: polyhedron support in VTK mesh scattering

ffed9d2

bd713 force-pushed the feat/bd713/meshScatter branch from 5b8df17 to ffed9d2 Compare July 28, 2026 02:03

joshua-white approved these changes Jul 28, 2026

View reviewed changes

bd713 and others added 2 commits August 11, 2026 14:07

Merge branch 'develop' into feat/bd713/meshScatter

a43e840

Added warning

e68ede0

victorapm removed the ci: run CUDA builds label Aug 27, 2026

bd713 and others added 2 commits September 17, 2026 14:07

Merge branch 'develop' into feat/bd713/meshScatter

70b742c

# Conflicts: # src/coreComponents/integrationTests/meshTests/testVTKImport.cpp # src/coreComponents/mesh/generators/VTKUtilities.cpp

Fixed vtk9.7

c3ae7bd

victorapm added the ci: run device builds Allows to run the device (CUDA/HIP) builds label Sep 18, 2026

Merge branch 'develop' into feat/bd713/meshScatter

deaa91b

joshua-white approved these changes Oct 22, 2026

View reviewed changes

bd713 and others added 6 commits October 22, 2026 16:47

Merge branch 'develop' into feat/bd713/meshScatter

e76bcf8

Fixed missing partitionRefinement

c7232fb

Merge remote-tracking branch 'origin/develop' into feat/bd713/meshSca…

5a6e2ae

…tter

schema

53643eb

Merge branch 'develop' into feat/bd713/meshScatter

ccd004a

Rebaseline

283936d

jhuang2601 approved these changes Oct 26, 2026

View reviewed changes

jafranc approved these changes Oct 28, 2026

View reviewed changes

-  for( vtkIdType i = 0; i < static_cast< vtkIdType >( assignment.size() ); ++i )
-  {
-    if( assignment[i] < mid )
-    {
-      loCells.push_back( i );
-      loAssign.push_back( assignment[i] );
-    }
-    else
-    {
-      hiCells.push_back( i );
-      hiAssign.push_back( assignment[i] );
-    }
-  }
+auto const n = static_cast< vtkIdType >( assignment.size() );
+stdVector< vtkIdType > order( n );
+std::iota( order.begin(), order.end(), vtkIdType{ 0 } );
+auto const split = std::stable_partition( order.begin(), order.end(),
+                                          [&]( vtkIdType i ) { return assignment[i] < mid; } );
+loCells.assign( order.begin(), split );
+hiCells.assign( split, order.end() );
+auto const rankOf = [&]( vtkIdType i ) { return assignment[i]; };
+std::ranges::transform( loCells, std::back_inserter( loAssign ), rankOf );
+std::ranges::transform( hiCells, std::back_inserter( hiAssign ), rankOf );

		std::function< void( vtkIdType, vtkIdType, integer, integer ) > bisect;
		bisect = [&]( vtkIdType begin, vtkIdType end, integer rankLo, integer rankHi )

	std::function< void( vtkIdType, vtkIdType, integer, integer ) > bisect;
	bisect = [&]( vtkIdType begin, vtkIdType end, integer rankLo, integer rankHi )
	auto bisect = [&]( vtkIdType begin, vtkIdType end, integer rankLo, integer rankHi )

Conversation

bd713 commented Apr 21, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Problem

Solution

Available strategies

Remark

Benchmark

Rebaseline needed

Uh oh!

victorapm commented Apr 23, 2026

Uh oh!

dkachuma left a comment

Choose a reason for hiding this comment

Uh oh!

joshua-white left a comment

Choose a reason for hiding this comment

Uh oh!

victorapm commented Oct 27, 2026

Uh oh!

jafranc left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Reviewers

Uh oh!

Assignees

Labels

Projects

Milestone

Development

Uh oh!

7 participants

bd713 commented Apr 21, 2026 •

edited

Loading