From 09b55fedee1e33c6054d398a8f1167be447ce3cc Mon Sep 17 00:00:00 2001
From: Soeren Peters <peters@irmb.tu-bs.de>
Date: Tue, 29 May 2018 16:36:01 +0200
Subject: [PATCH] - bugfix: find last or first fluid node to set BC

---
 .../grid/BoundaryConditions/Side.cpp          | 10 ++++--
 src/GridGenerator/grid/Field.cu               |  4 +--
 src/GridGenerator/grid/Grid.h                 |  3 +-
 .../grid/GridBuilder/LevelGridBuilder.cpp     |  7 +++-
 .../grid/GridBuilder/MultipleGridBuilder.cpp  | 34 +++++++++----------
 src/GridGenerator/grid/GridImp.cu             | 13 +++++++
 src/GridGenerator/grid/GridImp.h              |  6 ++--
 src/GridGenerator/grid/GridMocks.h            |  2 ++
 targets/apps/HULC/main.cpp                    |  9 +++--
 9 files changed, 60 insertions(+), 28 deletions(-)

diff --git a/src/GridGenerator/grid/BoundaryConditions/Side.cpp b/src/GridGenerator/grid/BoundaryConditions/Side.cpp
index d69cf285f..ec9053f56 100644
--- a/src/GridGenerator/grid/BoundaryConditions/Side.cpp
+++ b/src/GridGenerator/grid/BoundaryConditions/Side.cpp
@@ -105,7 +105,10 @@ void MX::addIndices(SPtr<Grid> grid, SPtr<BoundaryCondition> boundaryCondition)
     real startOuter = grid->getStartZ();
     real endOuter = grid->getEndZ();
 
-    Side::addIndices(grid, boundaryCondition, "x", grid->getStartX() + grid->getDelta(), startInner,
+    real coords[3] = { grid->getStartX(), grid->getStartY() + (grid->getEndY() - grid->getStartY()) / 2.0, grid->getStartZ() + (grid->getEndZ() - grid->getStartZ()) / 2.0 };
+    real startCoord = grid->getFirstFluidNode(coords, 0, grid->getStartX());
+
+    Side::addIndices(grid, boundaryCondition, "x", startCoord, startInner,
         endInner, startOuter, endOuter);
 }
 
@@ -119,7 +122,10 @@ void PX::addIndices(SPtr<Grid> grid, SPtr<BoundaryCondition> boundaryCondition)
     real startOuter = grid->getStartZ();
     real endOuter = grid->getEndZ();
 
-    Side::addIndices(grid, boundaryCondition, "x", grid->getEndX() - grid->getDelta(), startInner,
+    real coords[3] = { grid->getEndX(), grid->getStartY() + (grid->getEndY() - grid->getStartY()) / 2.0, grid->getStartZ() + (grid->getEndZ() - grid->getStartZ()) / 2.0 };
+    real startCoord = grid->getLastFluidNode(coords, 0, grid->getEndX());
+
+    Side::addIndices(grid, boundaryCondition, "x", startCoord, startInner,
         endInner, startOuter, endOuter);
 }
 
diff --git a/src/GridGenerator/grid/Field.cu b/src/GridGenerator/grid/Field.cu
index d03b51954..05949cef4 100644
--- a/src/GridGenerator/grid/Field.cu
+++ b/src/GridGenerator/grid/Field.cu
@@ -63,7 +63,7 @@ HOSTDEVICE bool Field::isFineToCoarseNode(uint index) const
 HOSTDEVICE bool Field::isFluid(uint index) const
 {
     const char type = field[index];
-    return type == FLUID || type == FLUID_CFC || type == FLUID_CFF || type == FLUID_FCC || type == FLUID_FCF || type == Q;
+    return type == FLUID || type == FLUID_CFC || type == FLUID_CFF || type == FLUID_FCC || type == FLUID_FCF || isRb(index);
 }
 
 HOSTDEVICE bool Field::isSolid(uint index) const
@@ -103,7 +103,7 @@ HOSTDEVICE bool Field::isQ(uint index) const
 
 HOSTDEVICE bool Field::isRb(uint index) const
 {
-    return  field[index] == SOLID;
+    return  field[index] == BC_GEOMETRY || field[index] == BC_OUTFLOW || field[index] == BC_VELOCITY || field[index] == BC_PRESSURE || field[index] == BC_SLIP;
 }
 
 // --------------------------------------------------------- //
diff --git a/src/GridGenerator/grid/Grid.h b/src/GridGenerator/grid/Grid.h
index cd71e64e6..7ac05f911 100644
--- a/src/GridGenerator/grid/Grid.h
+++ b/src/GridGenerator/grid/Grid.h
@@ -85,7 +85,8 @@ public:
 
     HOST virtual void findSparseIndices(SPtr<Grid> fineGrid) = 0;
 
-
+    HOSTDEVICE virtual real getFirstFluidNode(real coords[3], int direction, real startCoord) const = 0;
+    HOSTDEVICE virtual real getLastFluidNode(real coords[3], int direction, real startCoord) const = 0;
 };
 
 #endif
diff --git a/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.cpp b/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.cpp
index b0a3fb184..d8d164ed9 100644
--- a/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.cpp
+++ b/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.cpp
@@ -62,14 +62,15 @@ void LevelGridBuilder::setVelocityBoundaryCondition(SideType sideType, real vx,
         setVelocityGeometryBoundaryCondition(vx, vy, vz);
     else
     {
-
         SPtr<VelocityBoundaryCondition> velocityBoundaryCondition = VelocityBoundaryCondition::make(vx, vy, vz);
 
         auto side = SideFactory::make(sideType);
 
         velocityBoundaryConditions.push_back(velocityBoundaryCondition);
         velocityBoundaryCondition->side = side;
+        velocityBoundaryCondition->side->addIndices(grids[0], velocityBoundaryCondition);
     }
+
 }
 
 void LevelGridBuilder::setVelocityGeometryBoundaryCondition(real vx, real vy, real vz)
@@ -78,6 +79,7 @@ void LevelGridBuilder::setVelocityGeometryBoundaryCondition(real vx, real vy, re
     geometryBoundaryCondition->vx = vx;
     geometryBoundaryCondition->vy = vy;
     geometryBoundaryCondition->vz = vz;
+    geometryBoundaryCondition->side->addIndices(grids[0], geometryBoundaryCondition);
 }
 
 void LevelGridBuilder::setPressureBoundaryCondition(SideType sideType, real rho)
@@ -88,6 +90,7 @@ void LevelGridBuilder::setPressureBoundaryCondition(SideType sideType, real rho)
 
     pressureBoundaryConditions.push_back(pressureBoundaryCondition);
     pressureBoundaryCondition->side = side;
+    pressureBoundaryCondition->side->addIndices(grids[0], pressureBoundaryCondition);
 }
 
 void LevelGridBuilder::setPeriodicBoundaryCondition(bool periodic_X, bool periodic_Y, bool periodic_Z)
@@ -103,6 +106,8 @@ void LevelGridBuilder::setNoSlipBoundaryCondition(SideType sideType)
 
     noSlipBoundaryConditions.push_back(noSlipBoundaryCondition);
     noSlipBoundaryCondition->side = side;
+    noSlipBoundaryCondition->side->addIndices(grids[0], noSlipBoundaryCondition);
+
 }
 
 
diff --git a/src/GridGenerator/grid/GridBuilder/MultipleGridBuilder.cpp b/src/GridGenerator/grid/GridBuilder/MultipleGridBuilder.cpp
index 0034b71fe..51fefe703 100644
--- a/src/GridGenerator/grid/GridBuilder/MultipleGridBuilder.cpp
+++ b/src/GridGenerator/grid/GridBuilder/MultipleGridBuilder.cpp
@@ -271,23 +271,23 @@ void MultipleGridBuilder::buildGrids()
 
 
 
-    for(auto velocityBC : velocityBoundaryConditions)
-    {
-        velocityBC->side->addIndices(grids[0], velocityBC);
-    }
-
-    for (auto pressureBC : pressureBoundaryConditions)
-    {
-        pressureBC->side->addIndices(grids[0], pressureBC);
-    }
-
-    for (auto noSlipBC : noSlipBoundaryConditions)
-    {
-        noSlipBC->side->addIndices(grids[0], noSlipBC);
-    }
-
-    if(geometryBoundaryCondition)
-        geometryBoundaryCondition->side->addIndices(grids[0], geometryBoundaryCondition);
+    //for(auto velocityBC : velocityBoundaryConditions)
+    //{
+    //    velocityBC->side->addIndices(grids[0], velocityBC);
+    //}
+
+    //for (auto pressureBC : pressureBoundaryConditions)
+    //{
+    //    pressureBC->side->addIndices(grids[0], pressureBC);
+    //}
+
+    //for (auto noSlipBC : noSlipBoundaryConditions)
+    //{
+    //    noSlipBC->side->addIndices(grids[0], noSlipBC);
+    //}
+
+    //if(geometryBoundaryCondition)
+    //    geometryBoundaryCondition->side->addIndices(grids[0], geometryBoundaryCondition);
 
 }
 
diff --git a/src/GridGenerator/grid/GridImp.cu b/src/GridGenerator/grid/GridImp.cu
index 43f600853..c44b21000 100644
--- a/src/GridGenerator/grid/GridImp.cu
+++ b/src/GridGenerator/grid/GridImp.cu
@@ -501,6 +501,19 @@ HOSTDEVICE real GridImp::getNeighborCoord(bool periodicity, real startCoord, rea
     return coords[direction] + delta;
 }
 
+
+HOSTDEVICE real GridImp::getLastFluidNode(real coords[3], int direction, real startCoord) const
+{
+    coords[direction] = startCoord;
+    int index = this->transCoordToIndex(coords[0], coords[1], coords[2]);
+    while (!field.isFluid(index))
+    {
+        coords[direction] -= delta;
+        index = this->transCoordToIndex(coords[0], coords[1], coords[2]);
+    }
+    return coords[direction];
+}
+
 HOSTDEVICE real GridImp::getFirstFluidNode(real coords[3], int direction, real startCoord) const
 {
     coords[direction] = startCoord;
diff --git a/src/GridGenerator/grid/GridImp.h b/src/GridGenerator/grid/GridImp.h
index dee4aeccc..411886e3b 100644
--- a/src/GridGenerator/grid/GridImp.h
+++ b/src/GridGenerator/grid/GridImp.h
@@ -159,12 +159,14 @@ public:
 
     HOST void updateSparseIndices();
     HOSTDEVICE void setNeighborIndices(uint index);
-
+    HOSTDEVICE real getFirstFluidNode(real coords[3], int direction, real startCoord) const;
+    HOSTDEVICE real getLastFluidNode(real coords[3], int direction, real startCoord) const;
 private:
     HOSTDEVICE void setStopperNeighborCoords(uint index);
     HOSTDEVICE void getNeighborCoords(real &neighborX, real &neighborY, real &neighborZ, real x, real y, real z) const;
     HOSTDEVICE real getNeighborCoord(bool periodicity, real endCoord, real coords[3], int direction) const;
-    HOSTDEVICE real getFirstFluidNode(real coords[3], int direction, real startCoord) const;
+    
+
     HOSTDEVICE int getSparseIndex(const real &expectedX, const real &expectedY, const real &expectedZ) const;
 
     HOSTDEVICE static real getMinimumOnNodes(const real& minExact, const real& decimalStart, const real& delta);
diff --git a/src/GridGenerator/grid/GridMocks.h b/src/GridGenerator/grid/GridMocks.h
index 07e1cc6c6..0229958ad 100644
--- a/src/GridGenerator/grid/GridMocks.h
+++ b/src/GridGenerator/grid/GridMocks.h
@@ -72,6 +72,8 @@ public:
     void setPeriodicityZ(bool periodicity) override {}
     void findQs(Object* object) override {}
     void setFieldEntry(uint matrixIndex, char type) override {}
+    real getFirstFluidNode(real coords[3], int direction, real startCoord) const override { return 0.0; }
+    real getLastFluidNode(real coords[3], int direction, real startCoord) const override { return 0.0; }
 };
 
 class GridStub : public GridDummy
diff --git a/targets/apps/HULC/main.cpp b/targets/apps/HULC/main.cpp
index 3ca046c62..5cc779f5e 100644
--- a/targets/apps/HULC/main.cpp
+++ b/targets/apps/HULC/main.cpp
@@ -281,6 +281,10 @@ void multipleLevel(const std::string& configPath)
 
     gridBuilder->addGeometry(triangularMesh);
 
+
+    gridBuilder->buildGrids(); // buildGrids() has to be called before setting the BCs!!!!
+
+
     gridBuilder->setVelocityBoundaryCondition(SideType::MX, 0.001, 0.0, 0.0);
     gridBuilder->setPressureBoundaryCondition(SideType::PX, 0.001);
 
@@ -319,11 +323,10 @@ void multipleLevel(const std::string& configPath)
 
     //gridBuilder->writeGridToVTK("D:/GRIDGENERATION/gridTest_level_2", 2);
 
-    gridBuilder->buildGrids();
     //SimulationFileWriter::write("D:/GRIDGENERATION/files/", gridBuilder, FILEFORMAT::ASCII);
 
-    //gridBuilder->writeGridsToVtk("D:/GRIDGENERATION/");
-    //gridBuilder->writeArrows("D:/arrows");
+    gridBuilder->writeGridsToVtk("D:/GRIDGENERATION/");
+    gridBuilder->writeArrows("D:/arrows");
 
 
     SPtr<Parameter> para = Parameter::make();
-- 
GitLab