From f4d9d85102eee9c2cc4168536cf38ed5c9ec3bea Mon Sep 17 00:00:00 2001
From: Christoph Rettinger <christoph.rettinger@fau.de>
Date: Fri, 1 Dec 2017 10:11:45 +0100
Subject: [PATCH] Fixed mapping of infinite pe bodies

---
 src/pe_coupling/mapping/BodyBBMapping.cpp     | 23 +++++++++++++++----
 src/pe_coupling/mapping/BodyMapping.h         |  2 ++
 .../momentum_exchange_method/BodyMapping.h    |  2 ++
 3 files changed, 23 insertions(+), 4 deletions(-)

diff --git a/src/pe_coupling/mapping/BodyBBMapping.cpp b/src/pe_coupling/mapping/BodyBBMapping.cpp
index 393f87781..c3726a574 100644
--- a/src/pe_coupling/mapping/BodyBBMapping.cpp
+++ b/src/pe_coupling/mapping/BodyBBMapping.cpp
@@ -39,9 +39,26 @@ CellInterval getCellBB( const pe::ConstBodyID body, const IBlock & block, Struct
    if( body->isFinite() )
    {
       blockStorage.getCellBBFromAABB( cellBB, body->getAABB(), blockStorage.getLevel(block) );
-   } else
+   }
+   else
    {
-      blockStorage.getCellBBFromAABB( cellBB, body->getAABB().getIntersection( blockStorage.getDomain() ), blockStorage.getLevel(block) );
+      // if body is infinite (global), its AABB is also infinite which then requires special treatment
+
+      auto level = blockStorage.getLevel(block);
+      const real_t dx = blockStorage.dx(level);
+      const real_t dy = blockStorage.dy(level);
+      const real_t dz = blockStorage.dz(level);
+      Vector3<real_t> aabbExtensionByGhostLayers(real_c(numberOfGhostLayersToInclude) * dx,
+                                                 real_c(numberOfGhostLayersToInclude) * dy,
+                                                 real_c(numberOfGhostLayersToInclude) * dz);
+      auto extendedBlockAABB = blockStorage.getAABB(block.getId()).getExtended( aabbExtensionByGhostLayers );
+
+      // intersect the infinte (global) body with the block AABB, extended by its ghost layers
+      // then determine the cell bounding box of the intersection
+      blockStorage.getCellBBFromAABB( cellBB, body->getAABB().getIntersection( extendedBlockAABB ), level );
+
+      // if infinte body does not intersect with the extended block AABB, return an empty interval
+      if( cellBB.empty() ) return CellInterval();
    }
 
    cellBB.xMin() -= cell_idx_t(1); cellBB.yMin() -= cell_idx_t(1); cellBB.zMin() -= cell_idx_t(1);
@@ -61,7 +78,5 @@ CellInterval getCellBB( const pe::ConstBodyID body, const IBlock & block, Struct
    return cellBB;
 }
 
-
-
 } // namespace pe_coupling
 } // namespace walberla
diff --git a/src/pe_coupling/mapping/BodyMapping.h b/src/pe_coupling/mapping/BodyMapping.h
index 3c689ebb2..2e70c1439 100644
--- a/src/pe_coupling/mapping/BodyMapping.h
+++ b/src/pe_coupling/mapping/BodyMapping.h
@@ -54,6 +54,8 @@ void mapBody( const pe::BodyID & body, IBlock & block, StructuredBlockStorage &
 
    CellInterval cellBB = getCellBB( body, block, blockStorage, flagField->nrOfGhostLayers() );
 
+   if( cellBB.empty() ) return;
+
    Vector3<real_t> startCellCenter = blockStorage.getBlockLocalCellCenter( block, cellBB.min() );
    const real_t dx = blockStorage.dx( blockStorage.getLevel(block) );
    const real_t dy = blockStorage.dy( blockStorage.getLevel(block) );
diff --git a/src/pe_coupling/momentum_exchange_method/BodyMapping.h b/src/pe_coupling/momentum_exchange_method/BodyMapping.h
index 6149ec635..77d5cf401 100644
--- a/src/pe_coupling/momentum_exchange_method/BodyMapping.h
+++ b/src/pe_coupling/momentum_exchange_method/BodyMapping.h
@@ -198,6 +198,8 @@ void mapMovingBody( const pe::BodyID body, IBlock & block, StructuredBlockStorag
 
    CellInterval cellBB = getCellBB( body, block, blockStorage, flagField->nrOfGhostLayers() );
 
+   if( cellBB.empty() ) return;
+
    Vector3<real_t> startCellCenter = blockStorage.getBlockLocalCellCenter( block, cellBB.min() );
    const real_t dx = blockStorage.dx( blockStorage.getLevel(block) );
    const real_t dy = blockStorage.dy( blockStorage.getLevel(block) );
-- 
GitLab