reduce complexity of neighbor setup in migratePrimitives

before this was O(#neighbour_id * (#local_id + log(#neighbor_id)), now it should be O(#local_id + #neighbor_id)

src/hyteg/primitivestorage/PrimitiveStorage.cpp

@@ -41,6 +41,23 @@
 #include "hyteg/primitivestorage/SetupPrimitiveStorage.hpp"
 #include "hyteg/primitivestorage/loadbalancing/SimpleBalancer.hpp"
 
+
+namespace std{
+template<> struct hash<hyteg::PrimitiveID>{
+   size_t operator()(const hyteg::PrimitiveID &pid) const{
+      return std::hash<hyteg::PrimitiveID::IDType>{}(pid.getID());
+   }
+};
+
+template<> struct equal_to<hyteg::PrimitiveID>{
+   bool operator()( const hyteg::PrimitiveID& lhs, const hyteg::PrimitiveID& rhs ) const
+   {
+      return lhs.getID() == rhs.getID();
+   }
+};
+}
+
+
 namespace hyteg {
 
 using walberla::uint_t;
@@ -1136,67 +1153,69 @@ void PrimitiveStorage::migratePrimitives( const MigrationInfo& migrationInfo )
 
    std::vector< PrimitiveID > localPrimitiveIDs;
    getPrimitiveIDs( localPrimitiveIDs );
-   for ( const auto& localID : localPrimitiveIDs )
+
+   if(localPrimitiveIDs.empty())
    {
-      if ( vertexExistsInNeighborhood( localID ) )
-         neighborVertices_.erase( localID.getID() );
-      if ( edgeExistsInNeighborhood( localID ) )
-         neighborEdges_.erase( localID.getID() );
-      if ( faceExistsInNeighborhood( localID ) )
-         neighborFaces_.erase( localID.getID() );
-      if ( cellExistsInNeighborhood( localID ) )
-         neighborCells_.erase( localID.getID() );
+      neighborVertices_ = VertexMap{};
+      neighborEdges_ = EdgeMap{};
+      neighborFaces_ = FaceMap{};
+      neighborCells_ = CellMap{};
    }
+   else {
+      for ( const auto& localID : localPrimitiveIDs )
+      {
+         if ( vertexExistsInNeighborhood( localID ) )
+            neighborVertices_.erase( localID.getID() );
+         if ( edgeExistsInNeighborhood( localID ) )
+            neighborEdges_.erase( localID.getID() );
+         if ( faceExistsInNeighborhood( localID ) )
+            neighborFaces_.erase( localID.getID() );
+         if ( cellExistsInNeighborhood( localID ) )
+            neighborCells_.erase( localID.getID() );
+      }
 
-   /////////////////////////////////////////////////////////////////////////////////////////
-   // Erasing all neighbors that are not referenced by local primitives from neighborhood //
-   /////////////////////////////////////////////////////////////////////////////////////////
-
-   // enjoy worst possible complexity...
+      /////////////////////////////////////////////////////////////////////////////////////////
+      // Erasing all neighbors that are not referenced by local primitives from neighborhood //
+      /////////////////////////////////////////////////////////////////////////////////////////
 
-   std::vector< PrimitiveID > neighborhoodIDs;
-   for ( const auto& it : neighborVertices_ )
-      neighborhoodIDs.push_back( it.first );
-   for ( const auto& it : neighborEdges_ )
-      neighborhoodIDs.push_back( it.first );
-   for ( const auto& it : neighborFaces_ )
-      neighborhoodIDs.push_back( it.first );
-   for ( const auto& it : neighborCells_ )
-      neighborhoodIDs.push_back( it.first );
-
-   for ( const auto& neighborhoodID : neighborhoodIDs )
-   {
-      bool referenced = false;
+      // enjoy worst possible complexity...
+      std::vector< PrimitiveID > neighborhoodIDs;
+      neighborhoodIDs.reserve( neighborVertices_.size() + neighborEdges_.size() + neighborFaces_.size() + neighborCells_.size() );
+      for ( const auto& it : neighborVertices_ )
+         neighborhoodIDs.emplace_back( it.first );
+      for ( const auto& it : neighborEdges_ )
+         neighborhoodIDs.emplace_back( it.first );
+      for ( const auto& it : neighborFaces_ )
+         neighborhoodIDs.emplace_back( it.first );
+      for ( const auto& it : neighborCells_ )
+         neighborhoodIDs.emplace_back( it.first );
+
+      std::unordered_set<PrimitiveID> neighborhood_from_local_primitives;
+      std::vector< PrimitiveID > neighborIDs;
+      auto work = std::make_shared<std::vector< PrimitiveID >>();
       for ( const auto& localID : localPrimitiveIDs )
       {
-         if ( referenced )
-            break;
-
          Primitive* primitive = getPrimitive( localID );
 
-         std::vector< PrimitiveID > neighborIDs;
-         primitive->getNeighborPrimitives( neighborIDs );
+         neighborIDs.clear();
+         primitive->getNeighborPrimitives(neighborIDs, work);
 
-         for ( const auto& neighborOfLocalID : neighborIDs )
-         {
-            if ( neighborOfLocalID == neighborhoodID )
-            {
-               referenced = true;
-               break;
-            }
-         }
+         neighborhood_from_local_primitives.insert(neighborIDs.begin(), neighborIDs.end());
       }
 
-      if ( !referenced )
+      for ( const auto& neighborhoodID : neighborhoodIDs )
       {
-         if ( vertexExistsInNeighborhood( neighborhoodID ) )
-            neighborVertices_.erase( neighborhoodID.getID() );
-         if ( edgeExistsInNeighborhood( neighborhoodID ) )
-            neighborEdges_.erase( neighborhoodID.getID() );
-         if ( faceExistsInNeighborhood( neighborhoodID ) )
-            neighborFaces_.erase( neighborhoodID.getID() );
-         if ( cellExistsInNeighborhood( neighborhoodID ) )
-            neighborCells_.erase( neighborhoodID.getID() );
+         if ( neighborhood_from_local_primitives.find(neighborhoodID) == neighborhood_from_local_primitives.end())
+         {
+            if ( vertexExistsInNeighborhood( neighborhoodID ) )
+               neighborVertices_.erase( neighborhoodID.getID() );
+            if ( edgeExistsInNeighborhood( neighborhoodID ) )
+               neighborEdges_.erase( neighborhoodID.getID() );
+            if ( faceExistsInNeighborhood( neighborhoodID ) )
+               neighborFaces_.erase( neighborhoodID.getID() );
+            if ( cellExistsInNeighborhood( neighborhoodID ) )
+               neighborCells_.erase( neighborhoodID.getID() );
+         }
       }
    }