From ce57a46cf3deae0ce80541fc5e11b5549c6f7274 Mon Sep 17 00:00:00 2001
From: Soeren Peters <peters@irmb.tu-bs.de>
Date: Mon, 28 May 2018 15:45:02 +0200
Subject: [PATCH] - pressure BC add neighbor indices

---
 .../BoundaryConditions/BoundaryCondition.h    | 14 +++++++--
 .../grid/BoundaryConditions/Side.cpp          | 30 ++++++++++++++++++-
 .../grid/GridBuilder/GridBuilder.h            |  2 +-
 .../grid/GridBuilder/LevelGridBuilder.cpp     |  6 +++-
 .../grid/GridBuilder/LevelGridBuilder.h       |  2 +-
 .../GridReaderGenerator/GridGenerator.cpp     |  4 +--
 src/VirtualFluids_GPU/LBM/Simulation.cpp      | 12 ++++----
 targets/apps/HULC/main.cpp                    |  2 +-
 8 files changed, 57 insertions(+), 15 deletions(-)

diff --git a/src/GridGenerator/grid/BoundaryConditions/BoundaryCondition.h b/src/GridGenerator/grid/BoundaryConditions/BoundaryCondition.h
index 2250a4ea7..9b5e98f5e 100644
--- a/src/GridGenerator/grid/BoundaryConditions/BoundaryCondition.h
+++ b/src/GridGenerator/grid/BoundaryConditions/BoundaryCondition.h
@@ -12,6 +12,9 @@ class BoundaryCondition
 {
 public:
     std::vector<uint> indices;
+    SPtr<Side> side;
+
+    virtual void print() const = 0;
 
 };
 
@@ -23,13 +26,18 @@ public:
         return SPtr<PressureBoundaryCondition>(new PressureBoundaryCondition(rho));
     }
 
-    SPtr<Side> side;
+
+    std::vector<uint> neighborIndices;
+
     real rho;
 private:
     PressureBoundaryCondition(real rho) : rho(rho)
     {
 
     }
+
+public:
+    void print() const override {};
 };
 
 class VelocityBoundaryCondition : public BoundaryCondition
@@ -40,13 +48,15 @@ public:
         return SPtr<VelocityBoundaryCondition>(new VelocityBoundaryCondition(vx, vy, vz));
     }
 
-    SPtr<Side> side;
     real vx, vy, vz;
 private:
     VelocityBoundaryCondition(real vx, real vy, real vz) : vx(vx), vy(vy), vz(vz)
     {
 
     }
+
+public:
+    void print() const override {};
 };
 
 
diff --git a/src/GridGenerator/grid/BoundaryConditions/Side.cpp b/src/GridGenerator/grid/BoundaryConditions/Side.cpp
index e3c7d537c..ea7a73d36 100644
--- a/src/GridGenerator/grid/BoundaryConditions/Side.cpp
+++ b/src/GridGenerator/grid/BoundaryConditions/Side.cpp
@@ -12,13 +12,41 @@ void Side::addIndices(SPtr<Grid> grid, SPtr<BoundaryCondition> boundaryCondition
     {
         for (real v2 = startOuter; v2 <= endOuter; v2 += grid->getDelta())
         {
-            uint index = getIndex(grid, coord, constant, v1, v2);
+            const uint index = getIndex(grid, coord, constant, v1, v2);
             if (grid->getFieldEntry(index) == FLUID)
+            {
                 boundaryCondition->indices.push_back(grid->getSparseIndex(index) + 1);
+                setPressureNeighborIndices(boundaryCondition, grid, index);
+            }
         }
     }
 }
 
+void Side::setPressureNeighborIndices(SPtr<BoundaryCondition> boundaryCondition, SPtr<Grid> grid, const uint index)
+{
+    auto pressureBoundaryCondition = std::dynamic_pointer_cast<PressureBoundaryCondition>(boundaryCondition);
+    if (pressureBoundaryCondition)
+    {
+        real x, y, z;
+        grid->transIndexToCoords(index, x, y, z);
+
+        real nx = x;
+        real ny = y;
+        real nz = z;
+
+        if (boundaryCondition->side->getCoordinate() == X_INDEX)
+            nx = boundaryCondition->side->getDirection() * grid->getDelta() + x;
+        if (boundaryCondition->side->getCoordinate() == Y_INDEX)
+            ny = boundaryCondition->side->getDirection() * grid->getDelta() + y;
+        if (boundaryCondition->side->getCoordinate() == Z_INDEX)
+            nz = boundaryCondition->side->getDirection() * grid->getDelta() + z;
+
+        int neighborIndex = grid->transCoordToIndex(nx, ny, nz);
+        int sparseIndex = grid->getSparseIndex(neighborIndex);
+        pressureBoundaryCondition->neighborIndices.push_back(sparseIndex);
+    }
+}
+
 uint Side::getIndex(SPtr<Grid> grid, std::string coord, real constant, real v1, real v2)
 {
     if (coord == "x")
diff --git a/src/GridGenerator/grid/GridBuilder/GridBuilder.h b/src/GridGenerator/grid/GridBuilder/GridBuilder.h
index 1a5672c16..4f4a4e186 100644
--- a/src/GridGenerator/grid/GridBuilder/GridBuilder.h
+++ b/src/GridGenerator/grid/GridBuilder/GridBuilder.h
@@ -70,7 +70,7 @@ public:
     virtual uint getVelocitySize(int level) const = 0;
     virtual void getVelocityValues(real* vx, real* vy, real* vz, int* indices, int level) const = 0;
     virtual uint getPressureSize(int level) const = 0;
-    virtual void getPressureValues(real* rho, int* indices, int level) const = 0;
+    virtual void getPressureValues(real* rho, int* indices, int* neighborIndices, int level) const = 0;
     virtual void getVelocityQs(real* qs[27], int level) const = 0;
 };
 
diff --git a/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.cpp b/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.cpp
index 718f17561..a02f4c685 100644
--- a/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.cpp
+++ b/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.cpp
@@ -268,7 +268,7 @@ uint LevelGridBuilder::getPressureSize(int level) const
     return size;
 }
 
-void LevelGridBuilder::getPressureValues(real* rho, int* indices, int level) const
+void LevelGridBuilder::getPressureValues(real* rho, int* indices, int* neighborIndices, int level) const
 {
     int allIndicesCounter = 0;
     for (auto boundaryCondition : this->pressureBoundaryConditions)
@@ -276,12 +276,16 @@ void LevelGridBuilder::getPressureValues(real* rho, int* indices, int level) con
         for (int i = 0; i < boundaryCondition->indices.size(); i++)
         {
             indices[allIndicesCounter] = boundaryCondition->indices[i];
+
+            neighborIndices[allIndicesCounter] = boundaryCondition->neighborIndices[i];;
+
             rho[allIndicesCounter] = boundaryCondition->rho;
             allIndicesCounter++;
         }
     }
 }
 
+
 void LevelGridBuilder::getVelocityQs(real* qs[27], int level) const
 {
     int allIndicesCounter = 0;
diff --git a/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.h b/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.h
index 68f3636ed..add96e255 100644
--- a/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.h
+++ b/src/GridGenerator/grid/GridBuilder/LevelGridBuilder.h
@@ -70,7 +70,7 @@ public:
     VF_PUBLIC virtual void getVelocityValues(real* vx, real* vy, real* vz, int* indices, int level) const;
     VF_PUBLIC virtual void getVelocityQs(real* qs[27], int level) const;
     VF_PUBLIC uint getPressureSize(int level) const override;
-    VF_PUBLIC void getPressureValues(real* rho, int* indices, int level) const override;
+    VF_PUBLIC void getPressureValues(real* rho, int* indices, int* neighborIndices, int level) const override;
 
     VF_PUBLIC virtual void setPressValues(real* RhoBC, int* kN, int channelSide, int level) const;
 
diff --git a/src/VirtualFluids_GPU/DataStructureInitializer/GridReaderGenerator/GridGenerator.cpp b/src/VirtualFluids_GPU/DataStructureInitializer/GridReaderGenerator/GridGenerator.cpp
index 37d643992..7db737624 100644
--- a/src/VirtualFluids_GPU/DataStructureInitializer/GridReaderGenerator/GridGenerator.cpp
+++ b/src/VirtualFluids_GPU/DataStructureInitializer/GridReaderGenerator/GridGenerator.cpp
@@ -125,7 +125,7 @@ void GridGenerator::allocArrays_BoundaryValues()
 {
 	std::cout << "------read BoundaryValues------" << std::endl;
 
-    for (int level = 0; level < builder->getNumberOfGridLevels(); level++) {
+    for (uint level = 0; level < builder->getNumberOfGridLevels(); level++) {
         const auto numberOfPressureValues = int(builder->getPressureSize(level));
 
         cout << "size pressure level " << level << " : " << numberOfPressureValues << endl;
@@ -137,7 +137,7 @@ void GridGenerator::allocArrays_BoundaryValues()
         if (numberOfPressureValues > 1)
         {
             para->cudaAllocPress(level);
-            builder->getPressureValues(para->getParH(level)->QPress.RhoBC, para->getParH(level)->QPress.k, level);
+            builder->getPressureValues(para->getParH(level)->QPress.RhoBC, para->getParH(level)->QPress.k, para->getParH(level)->QPress.kN, level);
             para->cudaCopyPress(level);
         }
     }
diff --git a/src/VirtualFluids_GPU/LBM/Simulation.cpp b/src/VirtualFluids_GPU/LBM/Simulation.cpp
index 06e5521d1..928fd1f21 100644
--- a/src/VirtualFluids_GPU/LBM/Simulation.cpp
+++ b/src/VirtualFluids_GPU/LBM/Simulation.cpp
@@ -1586,12 +1586,12 @@ void Simulation::run()
 		  //getLastCudaError("QPressDevIncompNEQ27 execution failed");
 		  //////////////////////////////////////////////////////////////////////////////////
 		  //press NEQ comp
-		  //QPressDevNEQ27( para->getParD(0)->numberofthreads, para->getParD(0)->QPress.RhoBC, 
-		  //				  para->getParD(0)->d0SP.f[0],       para->getParD(0)->QPress.k,  
-		  //				  para->getParD(0)->QPress.kN,       para->getParD(0)->QPress.kQ,    para->getParD(0)->omega,
-		  //				  para->getParD(0)->neighborX_SP,    para->getParD(0)->neighborY_SP, para->getParD(0)->neighborZ_SP,
-		  //				  para->getParD(0)->size_Mat_SP,     para->getParD(0)->evenOrOdd);
-		  //getLastCudaError("QPressDevNEQ27 execution failed");
+          QPressDevNEQ27(para->getParD(0)->numberofthreads, para->getParD(0)->QPress.RhoBC,
+              para->getParD(0)->d0SP.f[0], para->getParD(0)->QPress.k,
+              para->getParD(0)->QPress.kN, para->getParD(0)->QPress.kQ, para->getParD(0)->omega,
+              para->getParD(0)->neighborX_SP, para->getParD(0)->neighborY_SP, para->getParD(0)->neighborZ_SP,
+              para->getParD(0)->size_Mat_SP, para->getParD(0)->evenOrOdd);
+          getLastCudaError("QPressDevNEQ27 execution failed");
 		  ////////////////////////////////////////////////////////////////////////////////
           //if (  myid == numprocs - 1)  
           //   PressSchlaffer27( para->getParD(0)->numberofthreads,  para->getParD(0)->Qoutflow.RhoBC,
diff --git a/targets/apps/HULC/main.cpp b/targets/apps/HULC/main.cpp
index e88f3aa4f..8022ad6a4 100644
--- a/targets/apps/HULC/main.cpp
+++ b/targets/apps/HULC/main.cpp
@@ -279,7 +279,7 @@ void multipleLevel(const std::string& configPath)
     //gridBuilder->addGeometry(triangularMesh);
 
     gridBuilder->setVelocityBoundaryCondition(SideType::MX, 0.001, 0.0, 0.0);
-    gridBuilder->setPressureBoundaryCondition(SideType::PX, 0.0);
+    gridBuilder->setPressureBoundaryCondition(SideType::PX, 0.001);
     //gridBuilder->setVelocityBoundaryCondition(SideType::PX, 0.001, 0.0, 0.0);
     //gridBuilder->setVelocityBoundaryCondition(SideType::MY, 0.001, 0.0, 0.0);
     //gridBuilder->setVelocityBoundaryCondition(SideType::PY, 0.001, 0.0, 0.0);
-- 
GitLab