From f19248f8849b253814712ce0a15cb4ddf797a466 Mon Sep 17 00:00:00 2001
From: Konstantin Kutscher <kutscher@irmb.tu-bs.de>
Date: Fri, 14 Jan 2022 14:40:44 +0100
Subject: [PATCH] fix undefine of variable interfaceWidth in LBMKernel, add
 variable phaseFieldBC for BC on the wall

---
 apps/cpu/LaminarTubeFlow/ltf.cfg              |  8 +--
 .../cpu/MultiphaseDropletTest/DropletTest.cfg |  2 +-
 apps/cpu/MultiphaseDropletTest/droplet.cpp    |  4 +-
 apps/cpu/RisingBubble2D/RisingBubble2D.cfg    | 12 ++--
 apps/cpu/RisingBubble2D/RisingBubble2D.cpp    | 61 ++++++++++---------
 .../CalculateTorqueCoProcessor.cpp            | 34 ++++++-----
 src/cpu/VirtualFluidsCore/LBM/LBMKernel.h     |  2 +-
 .../LBM/MultiphasePressureFilterLBMKernel.cpp |  9 +--
 .../LBM/MultiphasePressureFilterLBMKernel.h   | 12 ++++
 9 files changed, 77 insertions(+), 67 deletions(-)

diff --git a/apps/cpu/LaminarTubeFlow/ltf.cfg b/apps/cpu/LaminarTubeFlow/ltf.cfg
index 7ee23c7df..110b99d35 100644
--- a/apps/cpu/LaminarTubeFlow/ltf.cfg
+++ b/apps/cpu/LaminarTubeFlow/ltf.cfg
@@ -16,11 +16,11 @@ Re = 10
 
 logToFile = false
 
-newStart = false
+newStart = true
 restartStep = 10
 
-cpStart = 10
-cpStep = 10
+cpStart = 1000
+cpStep = 1000
 
 outTime = 10
-endTime = 10
\ No newline at end of file
+endTime = 1000
\ No newline at end of file
diff --git a/apps/cpu/MultiphaseDropletTest/DropletTest.cfg b/apps/cpu/MultiphaseDropletTest/DropletTest.cfg
index a93c46436..016e34072 100644
--- a/apps/cpu/MultiphaseDropletTest/DropletTest.cfg
+++ b/apps/cpu/MultiphaseDropletTest/DropletTest.cfg
@@ -1,5 +1,5 @@
 #pathname = d:/temp/MultiphaseDropletTest
-pathname = E:/Multiphase/DropletTest_new_corr5
+pathname = E:/Multiphase/DropletTest_Test
 
 numOfThreads = 4
 availMem = 10e9
diff --git a/apps/cpu/MultiphaseDropletTest/droplet.cpp b/apps/cpu/MultiphaseDropletTest/droplet.cpp
index a36636f23..49a887e98 100644
--- a/apps/cpu/MultiphaseDropletTest/droplet.cpp
+++ b/apps/cpu/MultiphaseDropletTest/droplet.cpp
@@ -73,7 +73,7 @@ void run(string configname)
 
 #if defined(__unix__)
          double lastTimeStep = 0;
-         if (!newStart) 
+         //if (!newStart) 
          {
              std::ifstream ifstr(fileName);
              ifstr >> lastTimeStep;
@@ -143,7 +143,7 @@ void run(string configname)
         //kernel = SPtr<LBMKernel>(new MultiphaseScratchCumulantLBMKernel());
        // kernel = SPtr<LBMKernel>(new MultiphaseCumulantLBMKernel());
         //kernel = SPtr<LBMKernel>(new MultiphaseTwoPhaseFieldsPressureFilterLBMKernel());
-        kernel = SPtr<LBMKernel>(new MultiphasePressureFilterCompressibleAirLBMKernel());
+        kernel = SPtr<LBMKernel>(new MultiphasePressureFilterLBMKernel());
 
         mu::Parser fgr;
         fgr.SetExpr("-(rho-rho_l)*g_y");
diff --git a/apps/cpu/RisingBubble2D/RisingBubble2D.cfg b/apps/cpu/RisingBubble2D/RisingBubble2D.cfg
index 99d68f3d0..659a3a2d2 100644
--- a/apps/cpu/RisingBubble2D/RisingBubble2D.cfg
+++ b/apps/cpu/RisingBubble2D/RisingBubble2D.cfg
@@ -1,4 +1,4 @@
-pathname = E:/Multiphase/RisingBubble2D_dr10
+pathname = E:/Multiphase/RisingBubble2D_dr10_test
 
 numOfThreads = 4
 availMem = 10e9
@@ -32,12 +32,12 @@ Mobility = 0.056 # 0.01 ./. 0.08, fine correction of Phase-field Relaxation para
 logToFile = false
 
 newStart = false
-restartStep = 50000
+restartStep = 10
 
-cpStart = 1000
-cpStep = 1000
+cpStart = 10
+cpStep = 10
 
-outTime = 100
-endTime = 100000
+outTime = 10000
+endTime = 20
 
 rStep = 159990 #160000
\ No newline at end of file
diff --git a/apps/cpu/RisingBubble2D/RisingBubble2D.cpp b/apps/cpu/RisingBubble2D/RisingBubble2D.cpp
index 939680185..e3a7e63ef 100644
--- a/apps/cpu/RisingBubble2D/RisingBubble2D.cpp
+++ b/apps/cpu/RisingBubble2D/RisingBubble2D.cpp
@@ -71,19 +71,19 @@ void run(string configname)
         
         std::string fileName = "./LastTimeStep" + std::to_string((int)boundingBox[1]) + ".txt";
 
-#if defined(__unix__)
-         double lastTimeStep = 0;
-         if (!newStart) 
-         {
-             std::ifstream ifstr(fileName);
-             ifstr >> lastTimeStep;
-             restartStep = lastTimeStep;
-             if(endTime >= lastTimeStep)
-                endTime = lastTimeStep + rStep;
-             else
-                return;
-         }    
-#endif
+//#if defined(__unix__)
+//         double lastTimeStep = 0;
+//         if (!newStart) 
+//         {
+//             std::ifstream ifstr(fileName);
+//             ifstr >> lastTimeStep;
+//             restartStep = lastTimeStep;
+//             if(endTime >= lastTimeStep)
+//                endTime = lastTimeStep + rStep;
+//             else
+//                return;
+//         }    
+//#endif
 
         //Sleep(30000);
 
@@ -165,6 +165,7 @@ void run(string configname)
         kernel->setDensityRatio(densityRatio);
         kernel->setMultiphaseModelParameters(beta, kappa);
         kernel->setContactAngle(theta);
+        dynamicPointerCast<MultiphasePressureFilterLBMKernel>(kernel)->setPhaseFieldBC(1.0);
 
         SPtr<BCProcessor> bcProc(new BCProcessor());
 
@@ -413,23 +414,23 @@ void run(string configname)
         if (myid == 0)
             UBLOG(logINFO, "Simulation-end");
             
-#if defined(__unix__)
-         //if (!newStart) 
-         //{
-            if (myid == 0) 
-            {
-                std::ofstream ostr(fileName);
-                ostr << endTime;
-                cout << "start sbatch\n";
-                //system("./start.sh");
-                //system("echo test!");
-                std::string str = "sbatch startJob" + std::to_string((int)boundingBox[1]) + ".sh";
-                //system("sbatch startJob512.sh");
-                system(str.c_str());
-            }   
-            //MPI_Barrier((MPI_Comm)comm->getNativeCommunicator()); 
-         //}
-#endif
+//#if defined(__unix__)
+//         //if (!newStart) 
+//         //{
+//            if (myid == 0) 
+//            {
+//                std::ofstream ostr(fileName);
+//                ostr << endTime;
+//                cout << "start sbatch\n";
+//                //system("./start.sh");
+//                //system("echo test!");
+//                std::string str = "sbatch startJob" + std::to_string((int)boundingBox[1]) + ".sh";
+//                //system("sbatch startJob512.sh");
+//                system(str.c_str());
+//            }   
+//            //MPI_Barrier((MPI_Comm)comm->getNativeCommunicator()); 
+//         //}
+//#endif
 
     } catch (std::exception &e) {
         cerr << e.what() << endl << flush;
diff --git a/src/cpu/VirtualFluidsCore/CoProcessors/CalculateTorqueCoProcessor.cpp b/src/cpu/VirtualFluidsCore/CoProcessors/CalculateTorqueCoProcessor.cpp
index fee4aaf3f..a3ccfd51f 100644
--- a/src/cpu/VirtualFluidsCore/CoProcessors/CalculateTorqueCoProcessor.cpp
+++ b/src/cpu/VirtualFluidsCore/CoProcessors/CalculateTorqueCoProcessor.cpp
@@ -103,16 +103,16 @@ void CalculateTorqueCoProcessor::calculateForces()
 
          SPtr<ILBMKernel> kernel = block->getKernel();
 
-         if (kernel->getCompressible())
-         {
-            calcMacrosFct = &D3Q27System::calcCompMacroscopicValues;
-            compressibleFactor = 1.0;
-         }
-         else
-         {
-            calcMacrosFct = &D3Q27System::calcIncompMacroscopicValues;
-            compressibleFactor = 0.0;
-         }
+         //if (kernel->getCompressible())
+         //{
+         //   calcMacrosFct = &D3Q27System::calcCompMacroscopicValues;
+         //   compressibleFactor = 1.0;
+         //}
+         //else
+         //{
+         //   calcMacrosFct = &D3Q27System::calcIncompMacroscopicValues;
+         //   compressibleFactor = 0.0;
+         //}
 
          SPtr<BCArray3D> bcArray = kernel->getBCProcessor()->getBCArray();          
          SPtr<DistributionArray3D> distributions = kernel->getDataSet()->getFdistributions(); 
@@ -160,10 +160,10 @@ void CalculateTorqueCoProcessor::calculateForces()
          //if we have got discretization with more level
          // deltaX is LBM deltaX and equal LBM deltaT 
          //double deltaX = 0.5; // LBMSystem::getDeltaT(block->getLevel()); //grid->getDeltaT(block);
-         double deltaXquadrat = deltaX*deltaX;
-         torqueX1 *= deltaXquadrat;
-         torqueX2 *= deltaXquadrat;
-         torqueX3 *= deltaXquadrat;
+         //double deltaXquadrat = deltaX*deltaX;
+         //torqueX1 *= deltaXquadrat;
+         //torqueX2 *= deltaXquadrat;
+         //torqueX3 *= deltaXquadrat;
 
          distributions->swap();
 
@@ -182,6 +182,8 @@ void CalculateTorqueCoProcessor::calculateForces()
    values.push_back(torqueX2global);
    values.push_back(torqueX3global);
 
+   //UBLOG(logINFO, "counter = " << counter);
+
    rvalues = comm->gather(values);
    if (comm->getProcessID() == comm->getRoot())
    {
@@ -204,8 +206,8 @@ UbTupleDouble3 CalculateTorqueCoProcessor::getForces(int x1, int x2, int x3,  SP
 
    LBMReal fs[D3Q27System::ENDF + 1];
    distributions->getDistributionInv(fs, x1, x2, x3);
-   LBMReal /*rho = 0.0,*/ vx1 = 0.0, vx2 = 0.0, vx3 = 0.0, drho = 0.0;
-   calcMacrosFct(fs, drho, vx1, vx2, vx3);
+   //LBMReal /*rho = 0.0,*/ vx1 = 0.0, vx2 = 0.0, vx3 = 0.0, drho = 0.0;
+   //calcMacrosFct(fs, drho, vx1, vx2, vx3);
    //rho = 1.0 + drho * compressibleFactor;
    
    if(bc)
diff --git a/src/cpu/VirtualFluidsCore/LBM/LBMKernel.h b/src/cpu/VirtualFluidsCore/LBM/LBMKernel.h
index 39205287b..1bdf610fd 100644
--- a/src/cpu/VirtualFluidsCore/LBM/LBMKernel.h
+++ b/src/cpu/VirtualFluidsCore/LBM/LBMKernel.h
@@ -165,7 +165,7 @@ protected:
     LBMReal phiH;
     LBMReal tauH;
     LBMReal mob;
-    LBMReal interfaceWidth;
+    LBMReal interfaceWidth { 4.0 };
 
 private:
     void checkFunction(mu::Parser fct);
diff --git a/src/cpu/VirtualFluidsCore/LBM/MultiphasePressureFilterLBMKernel.cpp b/src/cpu/VirtualFluidsCore/LBM/MultiphasePressureFilterLBMKernel.cpp
index f86e00896..0dc57aaba 100644
--- a/src/cpu/VirtualFluidsCore/LBM/MultiphasePressureFilterLBMKernel.cpp
+++ b/src/cpu/VirtualFluidsCore/LBM/MultiphasePressureFilterLBMKernel.cpp
@@ -87,6 +87,7 @@ SPtr<LBMKernel> MultiphasePressureFilterLBMKernel::clone()
 	kernel->setDeltaT(deltaT);
 	kernel->setGhostLayerWidth(2);
 	dynamicPointerCast<MultiphasePressureFilterLBMKernel>(kernel)->initForcing();
+    dynamicPointerCast<MultiphasePressureFilterLBMKernel>(kernel)->setPhaseFieldBC(this->phaseFieldBC);
 
 	return kernel;
 }
@@ -252,14 +253,8 @@ void MultiphasePressureFilterLBMKernel::calculate(int step)
 						+ (mfabb + mfcbb) + (mfbab + mfbcb) + (mfbba + mfbbc) + mfbbb;
 
 					LBMReal rho = rhoH + rhoToPhi * ((*phaseField)(x1, x2, x3) - phiH);
-					//! variable density -> TRANSFER!
-					//LBMReal rho = rhoH * ((*phaseField)(x1, x2, x3)) + rhoL * ((*phaseField2)(x1, x2, x3));
 
 					(*pressureOld)(x1, x2, x3) = (*pressure)(x1, x2, x3) + rho * c1o3 * drho;
-
-					//(*pressure)(x1, x2, x3) = (((*phaseField)(x1, x2, x3)) + ((*phaseField2)(x1, x2, x3)) - c1) * c1o3;
-					////!!!!!! relplace by pointer swap!
-					//(*pressureOld)(x1, x2, x3) = (*pressure)(x1, x2, x3);
 				}
 			}
 		}
@@ -1586,7 +1581,7 @@ void MultiphasePressureFilterLBMKernel::findNeighbors(CbArray3D<LBMReal, Indexer
 		if (!bcArray->isSolid(x1 + DX1[k], x2 + DX2[k], x3 + DX3[k])) {
 			phi[k] = (*ph)(x1 + DX1[k], x2 + DX2[k], x3 + DX3[k]);
 		} else {
-			phi[k] = 0.0;
+            phi[k] = phaseFieldBC;
 		}
 	}
 }
diff --git a/src/cpu/VirtualFluidsCore/LBM/MultiphasePressureFilterLBMKernel.h b/src/cpu/VirtualFluidsCore/LBM/MultiphasePressureFilterLBMKernel.h
index 05c4fd192..9b2b568b2 100644
--- a/src/cpu/VirtualFluidsCore/LBM/MultiphasePressureFilterLBMKernel.h
+++ b/src/cpu/VirtualFluidsCore/LBM/MultiphasePressureFilterLBMKernel.h
@@ -52,6 +52,16 @@ public:
     void calculate(int step) override;
     SPtr<LBMKernel> clone() override;
     double getCalculationTime() override { return .0; }
+
+    void setPhaseFieldBC(LBMReal bc)
+    {
+        phaseFieldBC = bc;
+    }
+    LBMReal getPhaseFieldBC()
+    {
+        return phaseFieldBC;
+    }
+
 protected:
     virtual void initDataSet();
     void swapDistributions() override;
@@ -94,6 +104,8 @@ protected:
     LBMReal forcingX1;
     LBMReal forcingX2;
     LBMReal forcingX3;
+
+    LBMReal phaseFieldBC { 0.0 }; // if 0.0 then light fluid on the wall, else if 1.0 havy fluid
 };
 
 #endif
-- 
GitLab