From 2266acb92997a757c2f62e5c3901b130f18570fc Mon Sep 17 00:00:00 2001
From: Kutscher <kutscher@irmb.tu-bs.de>
Date: Wed, 12 Jul 2023 09:46:46 +0200
Subject: [PATCH] fix bc issue in shotcrete jet setup
---
apps/cpu/ShotcreteJet/jet.cpp | 112 +++++++++++-------
.../Parallel/MPIIODataStructures.h | 14 +++
.../MPIIOMigrationSimulationObserver.cpp | 18 ++-
3 files changed, 97 insertions(+), 47 deletions(-)
diff --git a/apps/cpu/ShotcreteJet/jet.cpp b/apps/cpu/ShotcreteJet/jet.cpp
index 08482c956..f06fdefcd 100644
--- a/apps/cpu/ShotcreteJet/jet.cpp
+++ b/apps/cpu/ShotcreteJet/jet.cpp
@@ -97,7 +97,7 @@ int main(int argc, char *argv[])
double g_minX1 = -1.49631;
double g_minX2 = 0.193582;
- double g_minX3 = g_maxX3_box - 0.39;//-0.095; //-0.215;
+ double g_minX3 = g_maxX3_box - 0.03;//-0.095; //-0.215;
double g_maxX1 = -1.10631;
double g_maxX2 = 0.583582;
@@ -616,7 +616,7 @@ int main(int argc, char *argv[])
rcp->setBCSet(bcProc);
//////////////////////////////////////////////////////////////////////////
-if (newStart) {
+//if (newStart) {
SPtr<GbObject3D> gridCube = make_shared<GbCuboid3D>(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3);
if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), outputPath + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance());
@@ -725,7 +725,7 @@ if (newStart) {
///////////////////////////////////////////////////////////
SPtr<GbTriFaceMesh3D> meshInflowPipe2 = std::make_shared<GbTriFaceMesh3D>();
if (myid == 0) UBLOG(logINFO, "Read geoFluidArea:start");
- meshInflowPipe2->readMeshFromSTLFileBinary(geoPath + "/LongTube.stl", true);
+ meshInflowPipe2->readMeshFromSTLFileBinary(geoPath + "/LongTube2.stl", true);
if (myid == 0) UBLOG(logINFO, "Read geoFluidArea:end");
if (myid == 0) GbSystem3D::writeGeoObject(meshInflowPipe2.get(), outputPath + "/geo/LongTube", WbWriterVtkXmlBinary::getInstance());
SPtr<Interactor3D> intrInflowPipe2 = std::make_shared<D3Q27TriFaceMeshInteractor>(meshInflowPipe2, grid, noSlipBC, Interactor3D::SOLID, Interactor3D::EDGES);
@@ -734,23 +734,23 @@ if (newStart) {
///////////////////////////////////////////////////////////
//outflows
//////////////////////////////////////////////////////////
- SPtr<GbObject3D> geoOutflow1 = make_shared<GbCuboid3D>(g_minX1, g_minX2, g_minX3-2.0*dx, g_maxX1, g_maxX2, g_minX3);
+ SPtr<GbObject3D> geoOutflow1 = make_shared<GbCuboid3D>(g_minX1 - 2.0 * dx, g_minX2 - 2.0 * dx, g_minX3 - 2.0 * dx, g_maxX1 + 2.0 * dx, g_maxX2 + 2.0 * dx, g_minX3);
if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow1.get(), outputPath + "/geo/geoOutflow1", WbWriterVtkXmlBinary::getInstance());
SPtr<D3Q27Interactor> intrOutflow1 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow1, grid, outflowBC, Interactor3D::SOLID));
- SPtr<GbObject3D> geoOutflow2 = make_shared<GbCuboid3D>(g_minX1 - 2.0 * dx, g_minX2, g_minX3 - 2.0 * dx, g_minX1, g_maxX2, g_maxX3);
+ SPtr<GbObject3D> geoOutflow2 = make_shared<GbCuboid3D>(g_minX1 - 2.0 * dx, g_minX2 - 2.0 * dx, g_minX3 - 2.0 * dx, g_minX1, g_maxX2 + 2.0 * dx, g_maxX3 + 2.0 * dx);
if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow2.get(), outputPath + "/geo/geoOutflow2", WbWriterVtkXmlBinary::getInstance());
SPtr<D3Q27Interactor> intrOutflow2 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow2, grid, outflowBC, Interactor3D::SOLID));
- SPtr<GbObject3D> geoOutflow3 = make_shared<GbCuboid3D>(g_maxX1, g_minX2, g_minX3 - 2.0 * dx, g_maxX1 + 2.0 * dx, g_maxX2, g_maxX3);
+ SPtr<GbObject3D> geoOutflow3 = make_shared<GbCuboid3D>(g_maxX1, g_minX2 - 2.0 * dx, g_minX3 - 2.0 * dx, g_maxX1 + 2.0 * dx, g_maxX2 + 2.0 * dx, g_maxX3 + 2.0 * dx);
if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow3.get(), outputPath + "/geo/geoOutflow3", WbWriterVtkXmlBinary::getInstance());
SPtr<D3Q27Interactor> intrOutflow3 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow3, grid, outflowBC, Interactor3D::SOLID));
- SPtr<GbObject3D> geoOutflow4 = make_shared<GbCuboid3D>(g_minX1, g_minX2 - 2.0 * dx, g_minX3 - 2.0 * dx, g_maxX1, g_minX2, g_maxX3);
+ SPtr<GbObject3D> geoOutflow4 = make_shared<GbCuboid3D>(g_minX1 - 2.0 * dx, g_minX2 - 2.0 * dx, g_minX3 - 2.0 * dx, g_maxX1 + 2.0 * dx, g_minX2, g_maxX3 + 2.0 * dx);
if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow4.get(), outputPath + "/geo/geoOutflow4", WbWriterVtkXmlBinary::getInstance());
SPtr<D3Q27Interactor> intrOutflow4 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow4, grid, outflowBC, Interactor3D::SOLID));
- SPtr<GbObject3D> geoOutflow5 = make_shared<GbCuboid3D>(g_minX1, g_maxX2, g_minX3 - 2.0 * dx, g_maxX1, g_maxX2 + 2.0 * dx, g_maxX3);
+ SPtr<GbObject3D> geoOutflow5 = make_shared<GbCuboid3D>(g_minX1 - 2.0 * dx, g_maxX2, g_minX3 - 2.0 * dx, g_maxX1 + 2.0 * dx, g_maxX2 + 2.0 * dx, g_maxX3 + 2.0 * dx);
if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow5.get(), outputPath + "/geo/geoOutflow5", WbWriterVtkXmlBinary::getInstance());
SPtr<D3Q27Interactor> intrOutflow5 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow5, grid, outflowBC, Interactor3D::SOLID));
@@ -876,6 +876,13 @@ if (newStart) {
}
}
+ if (!newStart)
+ {
+ rcp->readBlocks((int)restartStep);
+ grid->accept(metisVisitor);
+ rcp->readDataSet((int)restartStep);
+ grid->setTimeStep(restartStep);
+ }
InteractorsHelper intHelper2(grid, metisVisitor, false);
intHelper2.addInteractor(intrInflowPipe);
@@ -959,32 +966,6 @@ if (newStart) {
ppblocks->update(0);
ppblocks.reset();
- // InitDistributionsBlockVisitor initVisitor;
- // grid->accept(initVisitor);
-
- double x1c = -1.31431 + R;
- double x2c = 0.375582 + R;
- double Ri = 5;
- double x3c = 0.136 + Ri;
-
- //R = 0.2 - 0.145; // 0.078-0.04; // 0.2;
-
- mu::Parser fct1;
- // fct1.SetExpr(" 0.5 - 0.5 * tanh(2 * (sqrt((x1 - x1c) ^ 2 + (x2 - x2c) ^ 2 + (x3 - x3c) ^ 2) - radius) / interfaceThickness)");
- fct1.SetExpr(" 0.5 - 0.5 * tanh(2 * (sqrt((x1 - x1c) ^ 2 + (x2 - x2c) ^ 2 + (x3 - x3c) ^ 2) - radius) / interfaceThickness)");
- fct1.DefineConst("x1c", x1c);
- fct1.DefineConst("x2c", x2c);
- fct1.DefineConst("x3c", x3c);
- fct1.DefineConst("radius", Ri);
- fct1.DefineConst("interfaceThickness", interfaceThickness * dx);
-
- MultiphaseVelocityFormInitDistributionsBlockVisitor initVisitor;
- initVisitor.setPhi(fct1);
- //grid->accept(initVisitor);
-
- //InitDistributionsBlockVisitor initVisitor;
- grid->accept(initVisitor);
-
// boundary conditions grid
{
SPtr<UbScheduler> geoSch(new UbScheduler(1));
@@ -995,24 +976,63 @@ if (newStart) {
-}
-else
-{
- //double restartStep = 10;
- rcp->restart((int)restartStep);
- grid->setTimeStep(restartStep);
-
- if (myid == 0) UBLOG(logINFO, "Restart - end");
-}
+ if (newStart)
+ {
+ double x1c = -1.31431 + R;
+ double x2c = 0.375582 + R;
+ double Ri = 5;
+ double x3c = 0.136 + Ri;
+
+ mu::Parser fct1;
+ // fct1.SetExpr(" 0.5 - 0.5 * tanh(2 * (sqrt((x1 - x1c) ^ 2 + (x2 - x2c) ^ 2 + (x3 - x3c) ^ 2) - radius) / interfaceThickness)");
+ fct1.SetExpr(" 0.5 - 0.5 * tanh(2 * (sqrt((x1 - x1c) ^ 2 + (x2 - x2c) ^ 2 + (x3 - x3c) ^ 2) - radius) / interfaceThickness)");
+ fct1.DefineConst("x1c", x1c);
+ fct1.DefineConst("x2c", x2c);
+ fct1.DefineConst("x3c", x3c);
+ fct1.DefineConst("radius", Ri);
+ fct1.DefineConst("interfaceThickness", interfaceThickness * dx);
+
+ MultiphaseVelocityFormInitDistributionsBlockVisitor initVisitor;
+ initVisitor.setPhi(fct1);
+ grid->accept(initVisitor);
+ }
+ //else
+ //{
+ // //rcp->restart((int)restartStep);
+ // rcp->readBlocks((int)restartStep);
+ // grid->accept(metisVisitor);
+ // rcp->readDataSet((int)restartStep);
+ // grid->setTimeStep(restartStep);
+ //}
+
+
+
+
+//}
+//else
+//{
+// //double restartStep = 10;
+// rcp->restart((int)restartStep);
+// grid->setTimeStep(restartStep);
+//
+// //GbCylinder3DPtr geoInflow(new GbCylinder3D(-1.30181 + 0.0005, 0.390872 - 0.00229, g_maxX3 - 2.0 * dx, -1.30181 + 0.0005, 0.390872 - 0.00229, g_maxX3 + 2.0 * dx, 0.013));
+// //if (myid == 0) GbSystem3D::writeGeoObject(geoInflow.get(), outputPath + "/geo/geoInflow", WbWriterVtkXmlBinary::getInstance());
+// //SPtr<D3Q27Interactor> intrInflow = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow, grid, inflowConcreteBC, Interactor3D::SOLID));
+// //SetBcBlocksBlockVisitor v1(intrInflow);
+// //grid->accept(v1);
+// //intrInflow->initInteractor();
+//
+// if (myid == 0) UBLOG(logINFO, "Restart - end");
+//}
- grid->accept(bcVisitor);
+ grid->accept(bcVisitor);
//OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm);
TwoDistributionsDoubleGhostLayerSetConnectorsBlockVisitor setConnsVisitor(comm);
// ThreeDistributionsDoubleGhostLayerSetConnectorsBlockVisitor setConnsVisitor(comm);
grid->accept(setConnsVisitor);
- int numOfThreads = 1;
+ int numOfThreads = 18;
omp_set_num_threads(numOfThreads);
SPtr<UbScheduler> nupsSch = std::make_shared<UbScheduler>(10, 10, 100);
@@ -1023,7 +1043,7 @@ else
// SPtr<UbScheduler> visSch(new UbScheduler(1, 8700, 8800));
// visSch->addSchedule(1, 8700, 8800);
SPtr<WriteSharpInterfaceQuantitiesSimulationObserver> writeMQSimulationObserver(new WriteSharpInterfaceQuantitiesSimulationObserver(grid, visSch, outputPath, WbWriterVtkXmlBinary::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm));
- writeMQSimulationObserver->update(0);
+ //writeMQSimulationObserver->update(10);
//SPtr<WriteMacroscopicQuantitiesSimulationObserver> writeMQSimulationObserver(new WriteMacroscopicQuantitiesSimulationObserver(grid, visSch, outputPath, WbWriterVtkXmlBinary::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm));
//writeMQSimulationObserver->update(0);
diff --git a/src/cpu/VirtualFluidsCore/Parallel/MPIIODataStructures.h b/src/cpu/VirtualFluidsCore/Parallel/MPIIODataStructures.h
index e84d4b7cd..e0f526889 100644
--- a/src/cpu/VirtualFluidsCore/Parallel/MPIIODataStructures.h
+++ b/src/cpu/VirtualFluidsCore/Parallel/MPIIODataStructures.h
@@ -64,6 +64,13 @@ struct DataSetRestart {
real collFactorL; // for Multiphase model
real collFactorG; // for Multiphase model
real densityRatio;// for Multiphase model
+ real sigma; // for Multiphase model
+ real contactAngle; // for Multiphase model
+ real phiL; // for Multiphase model
+ real phiH; // for Multiphase model
+ real tauH; // for Multiphase model
+ real mob; // for Multiphase model
+ real interfaceWidth; // for Multiphase model
int x1;
int x2;
int x3;
@@ -82,6 +89,13 @@ struct DataSetMigration {
real collFactorL; // for Multiphase model
real collFactorG; // for Multiphase model
real densityRatio;// for Multiphase model
+ real sigma; // for Multiphase model
+ real contactAngle; // for Multiphase model
+ real phiL; // for Multiphase model
+ real phiH; // for Multiphase model
+ real tauH; // for Multiphase model
+ real mob; // for Multiphase model
+ real interfaceWidth; // for Multiphase model
int globalID;
int ghostLayerWidth;
bool compressible;
diff --git a/src/cpu/VirtualFluidsCore/SimulationObservers/MPIIOMigrationSimulationObserver.cpp b/src/cpu/VirtualFluidsCore/SimulationObservers/MPIIOMigrationSimulationObserver.cpp
index 860b3f02a..cb36362df 100644
--- a/src/cpu/VirtualFluidsCore/SimulationObservers/MPIIOMigrationSimulationObserver.cpp
+++ b/src/cpu/VirtualFluidsCore/SimulationObservers/MPIIOMigrationSimulationObserver.cpp
@@ -31,7 +31,7 @@ MPIIOMigrationSimulationObserver::MPIIOMigrationSimulationObserver(SPtr<Grid3D>
//------------------------- define MPI types ---------------------------------
MPI_Datatype typesDataSet[3] = { MPI_DOUBLE, MPI_INT, MPI_CHAR };
- int blocksDataSet[3] = { 5, 2, 2 };
+ int blocksDataSet[3] = { 12, 2, 2 };
MPI_Aint offsetsDatatSet[3], lbDataSet, extentDataSet;
offsetsDatatSet[0] = 0;
@@ -173,6 +173,14 @@ void MPIIOMigrationSimulationObserver::writeDataSet(int step)
dataSetArray[ic].collFactorG = kernel->getCollisionFactorG();
dataSetArray[ic].densityRatio = kernel->getDensityRatio();
+ dataSetArray[ic].sigma = kernel->getSigma();
+ dataSetArray[ic].contactAngle = kernel->getContactAngle();
+ dataSetArray[ic].phiL = kernel->getPhiL();
+ dataSetArray[ic].phiH = kernel->getPhiH();
+ dataSetArray[ic].tauH = kernel->getPhaseFieldRelaxation();
+ dataSetArray[ic].mob = kernel->getMobility();
+ dataSetArray[ic].interfaceWidth = kernel->getInterfaceWidth();
+
D3Q27EsoTwist3DSplittedVectorPtrF = dynamicPointerCast<D3Q27EsoTwist3DSplittedVector>(block->getKernel()->getDataSet()->getFdistributions());
localDistributionsF = D3Q27EsoTwist3DSplittedVectorPtrF->getLocalDistributions();
nonLocalDistributionsF = D3Q27EsoTwist3DSplittedVectorPtrF->getNonLocalDistributions();
@@ -1154,6 +1162,14 @@ void MPIIOMigrationSimulationObserver::readDataSet(int step)
kernel->setCollisionFactorMultiphase(dataSetArray[n].collFactorL, dataSetArray[n].collFactorG);
kernel->setDensityRatio(dataSetArray[n].densityRatio);
+ kernel->setSigma(dataSetArray[n].sigma);
+ kernel->setContactAngle(dataSetArray[n].contactAngle);
+ kernel->setPhiL(dataSetArray[n].phiL);
+ kernel->setPhiH(dataSetArray[n].phiH);
+ kernel->setPhaseFieldRelaxation(dataSetArray[n].tauH);
+ kernel->setMobility(dataSetArray[n].mob);
+ kernel->setInterfaceWidth(dataSetArray[n].interfaceWidth);
+
SPtr<DataSet3D> dataSetPtr = SPtr<DataSet3D>(new DataSet3D());
dataSetPtr->setFdistributions(mFdistributions);
if (multiPhase1)
--
GitLab