diff --git a/source/Applications/Thermoplast/config.txt b/source/Applications/Thermoplast/config.txt
index 392cff19fd6b2be6a89057b918a2ab7dfa0ad51b..7ddd18eac4382cfb26474f5a051d4c5266037949 100644
--- a/source/Applications/Thermoplast/config.txt
+++ b/source/Applications/Thermoplast/config.txt
@@ -1,21 +1,24 @@
#simulation parameters
#boundingBox = 0 0 0 300 1520 2320
boundingBox = 60 1370 10 190 1530 320
+#boundingBox = 60 1370 120 100 1530 320
+
+#boundingBox = 60 0 10 190 1530 2320
blocknx = 10 10 10
#blocknx = 300 420 320
-endTime = 100000
-outTime = 100
+endTime = 20
+outTime = 20
availMem = 25e9
-uLB = 0.01
+uLB = 0.1
#PE parameters
-peMinOffset = 50 2 2
-peMaxOffset = -8 -30 -2
-sphereTime = 1000
+peMinOffset = 0 0 0 #46 2 2
+peMaxOffset = 0 0 0 #-8 -25 -2
+sphereTime = 20
#geometry files
pathGeo = d:/Projects/ThermoPlast/SimGeo
michel = /michel.stl
plexiglas = /plexiglas.stl
-pathOut = d:/temp/thermoplastCluster
\ No newline at end of file
+pathOut = g:/temp/thermoplastCluster3
\ No newline at end of file
diff --git a/source/Applications/Thermoplast/thermoplast.cpp b/source/Applications/Thermoplast/thermoplast.cpp
index 24f26a6184e099791a61ab0b60a8ab89c899d8f8..275e5b8778e8ac8f5eada6fe2ee35a5fde55b498 100644
--- a/source/Applications/Thermoplast/thermoplast.cpp
+++ b/source/Applications/Thermoplast/thermoplast.cpp
@@ -31,6 +31,7 @@
#include <WriteDemObjectsCoProcessor.h>
#include "CreateDemObjectsCoProcessor.h"
+#include "RestartDemCouplingCoProcessor.h"
using namespace std;
@@ -46,10 +47,10 @@ double g_maxX3 = 0;
vector<double> peMinOffset;
vector<double> peMaxOffset;
-string pathOut = "d:/temp/thermoplastCluster";
-string pathGeo = "d:/Projects/ThermoPlast/Geometrie";
+string pathOut;// = "d:/temp/thermoplastCluster";
+string pathGeo;// = "d:/Projects/ThermoPlast/Geometrie";
-std::shared_ptr<DemCoProcessor> makePeCoProcessor(SPtr<Grid3D> grid, SPtr<Communicator> comm, const SPtr<UbScheduler> peScheduler, const std::shared_ptr<LBMUnitConverter> lbmUnitConverter, int maxpeIterations)
+std::shared_ptr<DemCoProcessor> makePeCoProcessor(SPtr<Grid3D> grid, SPtr<Communicator> comm, const SPtr<UbScheduler> peScheduler, const std::shared_ptr<LBMUnitConverter> lbmUnitConverter, int maxpeIterations)
{
double peRelaxtion = 0.7;
//int maxpeIterations = 10000;
@@ -65,17 +66,17 @@ std::shared_ptr<DemCoProcessor> makePeCoProcessor(SPtr<Grid3D> grid, SPtr<Commun
//UbTupleInt3 simulationDomain(gridNx, gridNy, gridNz);
//std::array<double, 6> simulationDomain = {1, 1, 1, 30, 30, 30};
- std::array<double, 6> simulationDomain = {g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3};
+ std::array<double, 6> simulationDomain ={ g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3 };
UbTupleInt3 numberOfBlocks(grid->getNX1(), grid->getNX2(), grid->getNX3());
//UbTupleInt3 numberOfBlocks((simulationDomain[3]-simulationDomain[0])/val<1>(grid->getBlockNX()), (simulationDomain[4]-simulationDomain[1])/val<2>(grid->getBlockNX()), (simulationDomain[5]-simulationDomain[2])/val<3>(grid->getBlockNX()));
UbTupleBool3 isPeriodic(grid->isPeriodicX1(), grid->isPeriodicX2(), grid->isPeriodicX3());
Vector3D minOffset(peMinOffset[0], peMinOffset[1], peMinOffset[2]);
Vector3D maxOffset(peMaxOffset[0], peMaxOffset[1], peMaxOffset[2]);
-
+
SPtr<GbObject3D> boxPE(new GbCuboid3D(simulationDomain[0]+minOffset[0], simulationDomain[1]+minOffset[1], simulationDomain[2]+minOffset[2], simulationDomain[3]+maxOffset[0], simulationDomain[4]+maxOffset[1], simulationDomain[5]+maxOffset[2]));
GbSystem3D::writeGeoObject(boxPE.get(), pathOut + "/geo/boxPE", WbWriterVtkXmlBinary::getInstance());
- std::shared_ptr<PeParameter> peParamter = std::make_shared<PeParameter>(peRelaxtion, maxpeIterations, globalLinearAcc,
+ std::shared_ptr<PeParameter> peParamter = std::make_shared<PeParameter>(peRelaxtion, maxpeIterations, globalLinearAcc,
planeMaterial, simulationDomain, numberOfBlocks, isPeriodic, minOffset, maxOffset);
std::shared_ptr<PhysicsEngineSolverAdapter> peSolver = std::make_shared<PePhysicsEngineSolverAdapter>(peParamter);
@@ -94,349 +95,387 @@ void thermoplast(string configname)
{
try
{
- SPtr<Communicator> comm = MPICommunicator::getInstance();
- int myid = comm->getProcessID();
-
- ConfigurationFile config;
- config.load(configname);
-
- vector<int> blocknx = config.getVector<int>("blocknx");
- vector<double> boundingBox = config.getVector<double>("boundingBox");
- peMinOffset = config.getVector<double>("peMinOffset");
- peMaxOffset = config.getVector<double>("peMaxOffset");
- int endTime = config.getValue<int>("endTime");
- double outTime = config.getValue<double>("outTime");
- double availMem = config.getValue<double>("availMem");
- double uLB = config.getValue<double>("uLB");
- string pathOut = config.getValue<string>("pathOut");
- string pathGeo = config.getValue<string>("pathGeo");
- string michel = config.getValue<string>("michel");
- string plexiglas = config.getValue<string>("plexiglas");
- double sphereTime = config.getValue<double>("sphereTime");
-
- //parameters
- //string pathOut = "d:/temp/thermoplast3";
- //string pathGeo = "d:/Projects/ThermoPlast/Geometrie";
- int numOfThreads = 1;
- //int blocknx[3] ={ 10,10,10 };
- //double endTime = 1000000;
- //double outTime = 300;
- //double availMem = 8e9;
- double deltax = 1;
- double rhoLB = 0.0;
- //double uLB = 0.1;
- double radius = 5;
- double Re = 300;
- double nuLB = (uLB*2.0*radius)/Re;
-
-
-
- //geometry definition
-
- //simulation bounding box
- g_minX1 = boundingBox[0];
- g_minX2 = boundingBox[1];
- g_minX3 = boundingBox[2];
-
- g_maxX1 = boundingBox[3];
- g_maxX2 = boundingBox[4];
- g_maxX3 = boundingBox[5];
-
- double blockLength = blocknx[0]*deltax;
-
- SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3));
- if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), pathOut + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance());
-
- //box
- SPtr<GbObject3D> box(new GbCuboid3D(g_minX1-blockLength, g_minX2, g_minX3, g_maxX1+blockLength, g_maxX2, g_maxX3));
- GbSystem3D::writeGeoObject(box.get(), pathOut + "/geo/box", WbWriterVtkXmlBinary::getInstance());
-
- //michel
- if (myid==0) UBLOG(logINFO, "Read michelGeo:start");
- SPtr<GbTriFaceMesh3D> michelGeo = SPtr<GbTriFaceMesh3D>(GbTriFaceMesh3DCreator::getInstance()->readMeshFromSTLFile2(pathGeo+michel, "michelGeo", GbTriFaceMesh3D::KDTREE_SAHPLIT, false));
- if (myid==0) UBLOG(logINFO, "Read michelGeo:end");
- if (myid==0) GbSystem3D::writeGeoObject(michelGeo.get(), pathOut+"/geo/michelGeo", WbWriterVtkXmlBinary::getInstance());
-
- //plexiglas
- if (myid==0) UBLOG(logINFO, "Read plexiglasGeo:start");
- SPtr<GbTriFaceMesh3D> plexiglasGeo = SPtr<GbTriFaceMesh3D>(GbTriFaceMesh3DCreator::getInstance()->readMeshFromSTLFile2(pathGeo+plexiglas, "plexiglasGeo", GbTriFaceMesh3D::KDTREE_SAHPLIT, false));
- if (myid==0) UBLOG(logINFO, "Read plexiglasGeo:end");
- if (myid==0) GbSystem3D::writeGeoObject(plexiglasGeo.get(), pathOut+"/geo/plexiglasGeo", WbWriterVtkXmlBinary::getInstance());
-
- if (myid == 0)
- {
- UBLOG(logINFO, "Parameters:");
- UBLOG(logINFO, "* uLB = " << uLB);
- UBLOG(logINFO, "* rhoLB = " << rhoLB);
- UBLOG(logINFO, "* nuLB = " << nuLB);
- UBLOG(logINFO, "* deltaX = " << deltax);
- UBLOG(logINFO, "* radius = " << radius);
- UBLOG(logINFO, "* Re = " << Re);
- UBLOG(logINFO, "* number of threads = "<<numOfThreads);
- UBLOG(logINFO, "* number of processes = "<<comm->getNumberOfProcesses());
- UBLOG(logINFO, "* path = "<<pathOut);
- UBLOG(logINFO, "Preprocess - start");
- }
-
- //Grid definition
- SPtr<Grid3D> grid(new Grid3D(comm));
- grid->setDeltaX(deltax);
- grid->setBlockNX(blocknx[0], blocknx[1], blocknx[2]);
- grid->setPeriodicX1(false);
- grid->setPeriodicX2(false);
- grid->setPeriodicX3(false);
-
- //blocks generating
- GenBlocksGridVisitor genBlocks(gridCube);
- grid->accept(genBlocks);
-
- //boundary conditions definition
- //boundary conditions adapters
- //////////////////////////////////////////////////////////////////////////////
- SPtr<BCAdapter> noSlipBCAdapter(new NoSlipBCAdapter());
- noSlipBCAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NoSlipBCAlgorithm()));
-
- mu::Parser fct;
- fct.SetExpr("U");
- fct.DefineConst("U", uLB);
- SPtr<BCAdapter> inflowAdapter(new VelocityBCAdapter(true, false, false, fct, 0, BCFunction::INFCONST));
- inflowAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityBCAlgorithm()));
-
- SPtr<BCAdapter> outflowAdapter(new DensityBCAdapter(rhoLB));
- outflowAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NonEqDensityBCAlgorithm()));
- //outflowAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NonReflectingOutflowBCAlgorithm()));
-
- //sphere
- mu::Parser fct2;
- fct2.SetExpr("U");
- fct2.DefineConst("U", 0.0);
- SPtr<BCAdapter> velocityBcParticleAdapter(new VelocityBCAdapter(true, false, false, fct2, 0, BCFunction::INFCONST));
- velocityBcParticleAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityWithDensityBCAlgorithm()));
-
- //inflow
- GbCuboid3DPtr geoInflow1(new GbCuboid3D(g_minX1-blockLength, g_maxX2-120.0, g_minX3+190.0, g_minX1+1, g_maxX2+20.0, g_minX3+320.0));
if (myid == 0) GbSystem3D::writeGeoObject(geoInflow1.get(), pathOut + "/geo/geoInflow1", WbWriterVtkXmlASCII::getInstance());
-
- GbCuboid3DPtr geoInflow2(new GbCuboid3D(g_minX1-blockLength, g_minX2 - blockLength, g_minX3 - blockLength, g_minX1, g_maxX2 + blockLength, g_maxX3 + blockLength));
- if (myid == 0) GbSystem3D::writeGeoObject(geoInflow2.get(), pathOut + "/geo/geoInflow2", WbWriterVtkXmlASCII::getInstance());
-
- GbCuboid3DPtr geoInflow3(new GbCuboid3D(g_minX1-blockLength, g_minX2-radius, g_maxX3-4.0*radius-1, g_minX1+1, g_maxX2+radius, g_maxX3-4.0*radius));
- if (myid == 0) GbSystem3D::writeGeoObject(geoInflow3.get(), pathOut + "/geo/geoInflow3", WbWriterVtkXmlASCII::getInstance());
-
- GbCuboid3DPtr geoInflow4(new GbCuboid3D(g_minX1-blockLength, g_minX2+4.0*radius, g_maxX3-4.0*radius-1.0, g_minX1+1, g_minX2+4.0*radius+1.0, g_maxX3+radius));
- if (myid == 0) GbSystem3D::writeGeoObject(geoInflow4.get(), pathOut + "/geo/geoInflow4", WbWriterVtkXmlASCII::getInstance());
-
- //outflow
- GbCuboid3DPtr geoOutflow(new GbCuboid3D(g_maxX1, g_minX2 - blockLength, g_minX3 - blockLength, g_maxX1 + blockLength, g_maxX2 + blockLength, g_maxX3 + blockLength));
- if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow.get(), pathOut + "/geo/geoOutflow", WbWriterVtkXmlASCII::getInstance());
-
- //boundary conditions visitor
- SPtr<BoundaryConditionsBlockVisitor> bcVisitor(new BoundaryConditionsBlockVisitor());
- bcVisitor->addBC(noSlipBCAdapter);
- bcVisitor->addBC(inflowAdapter);
- bcVisitor->addBC(outflowAdapter);
- bcVisitor->addBC(velocityBcParticleAdapter);
- //////////////////////////////////////////////////////////////////////////////////
-
- //set boundary conditions for blocks and create process decomposition for MPI
- SPtr<D3Q27Interactor> boxInt(new D3Q27Interactor(box, grid, noSlipBCAdapter, Interactor3D::INVERSESOLID));
-
- //sphere bc object
- //SPtr<MovableObjectInteractor> sphereInt1;
- //const std::shared_ptr<Reconstructor> velocityReconstructor = std::make_shared<VelocityBcReconstructor>();
- //const std::shared_ptr<Reconstructor> extrapolationReconstructor = std::make_shared<ExtrapolationReconstructor>(velocityReconstructor);
- //sphereInt1 = SPtr<MovableObjectInteractor>(new MovableObjectInteractor(sphere1, grid, velocityBcParticleAdapter, Interactor3D::SOLID, extrapolationReconstructor, State::UNPIN));
-
- //inflow
- SPtr<D3Q27Interactor> inflowInt1 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow1, grid, inflowAdapter, Interactor3D::SOLID));
- SPtr<D3Q27Interactor> inflowInt2 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow2, grid, outflowAdapter, Interactor3D::SOLID));
-
- SPtr<D3Q27Interactor> inflowInt3 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow3, grid, noSlipBCAdapter, Interactor3D::SOLID));
- SPtr<D3Q27Interactor> inflowInt4 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow4, grid, noSlipBCAdapter, Interactor3D::SOLID));
-
- //outflow
- SPtr<D3Q27Interactor> outflowInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow, grid, outflowAdapter, Interactor3D::SOLID));
-
- //michel
- SPtr<Interactor3D> michelInt = SPtr<D3Q27TriFaceMeshInteractor>(new D3Q27TriFaceMeshInteractor(michelGeo, grid, noSlipBCAdapter, Interactor3D::SOLID));
-
- //plexiglas
- SPtr<Interactor3D> plexiglasInt = SPtr<D3Q27TriFaceMeshInteractor>(new D3Q27TriFaceMeshInteractor(plexiglasGeo, grid, noSlipBCAdapter, Interactor3D::SOLID));
-
- //PE
- double refLengthLb = radius*2.0;
- double refLengthWorld = refLengthLb * deltax;
- const std::shared_ptr<LBMUnitConverter> lbmUnitConverter = std::make_shared<LBMUnitConverter>(refLengthWorld, LBMUnitConverter::WORLD_MATERIAL::AIR_20C, refLengthLb);
- if (myid == 0) std::cout << lbmUnitConverter->toString() << std::endl;
- double rhoSphere = 915 * lbmUnitConverter->getFactorDensityWToLb(); // kg/m^3
- if (myid == 0) UBLOG(logINFO, "rhoSphere = "<<rhoSphere);
- SPtr<PhysicsEngineMaterialAdapter> sphereMaterial(new PePhysicsEngineMaterialAdapter("Polypropylen", rhoSphere, 0, 0.15, 0.1, 0.45, 0.5, 1, 0, 0));
- const int timestep = 2;
- const SPtr<UbScheduler> peScheduler(new UbScheduler(timestep));
- int maxpeIterations = 10;//endTime/2;
- SPtr<DemCoProcessor> demCoProcessor = makePeCoProcessor(grid, comm, peScheduler, lbmUnitConverter, maxpeIterations);
- //demCoProcessor->addInteractor(sphereInt1, sphereMaterial, Vector3D(0.0, 0.0, 0.0));
- //demCoProcessor->distributeIDs();
- demCoProcessor->setBlockVisitor(bcVisitor);
- //double g = 9.81 * lbmUnitConverter->getFactorAccWToLb();
- //double rhoFluid = lbmUnitConverter->getFactorDensityWToLb() * 830.0; // 1 // kg/m^3
- //////////////////////////////////////////////////////////////////////////
- SPtr<Grid3DVisitor> peVisitor(new PePartitioningGridVisitor(comm, demCoProcessor));
- //grid->accept(peVisitor);
-
- //SetSolidBlocksBlockVisitor solidBlocksVisitor1(michelInt);
- //grid->accept(solidBlocksVisitor1);
- //SetSolidBlocksBlockVisitor solidBlocksVisitor2(plexiglasInt);
- //grid->accept(solidBlocksVisitor2);
-
- //SetBcBlocksBlockVisitor bcBlocksVisitor1(michelInt);
- //grid->accept(bcBlocksVisitor1);
- //SetBcBlocksBlockVisitor bcBlocksVisitor2(plexiglasInt);
- //grid->accept(bcBlocksVisitor2);
-
- InteractorsHelper intHelper(grid, peVisitor);
- intHelper.addInteractor(boxInt);
- //intHelper.addInteractor(sphereInt1);
- //intHelper.addInteractor(inflowInt2);
- intHelper.addInteractor(michelInt);
- intHelper.addInteractor(plexiglasInt);
- intHelper.addInteractor(inflowInt1);
- intHelper.addInteractor(outflowInt);
- intHelper.addInteractor(inflowInt2);
- //intHelper.addInteractor(inflowInt3);
- //intHelper.addInteractor(inflowInt4);
-
- //intHelper.addInteractor(plexiglasInt);
-
- intHelper.selectBlocks();
-
- //write data for visualization of block grid
- SPtr<CoProcessor> ppblocks(new WriteBlocksCoProcessor(grid, SPtr<UbScheduler>(new UbScheduler(1)), pathOut, WbWriterVtkXmlBinary::getInstance(), comm));
- ppblocks->process(0);
- ppblocks.reset();
-
- unsigned long long numberOfBlocks = (unsigned long long)grid->getNumberOfBlocks();
- int ghostLayer = 3;
- unsigned long long numberOfNodesPerBlock = (unsigned long long)(blocknx[0])* (unsigned long long)(blocknx[1])* (unsigned long long)(blocknx[2]);
- unsigned long long numberOfNodes = numberOfBlocks * numberOfNodesPerBlock;
- unsigned long long numberOfNodesPerBlockWithGhostLayer = numberOfBlocks * (blocknx[0] + ghostLayer) * (blocknx[1] + ghostLayer) * (blocknx[2] + ghostLayer);
- double needMemAll = double(numberOfNodesPerBlockWithGhostLayer*(27 * sizeof(double) + sizeof(int) + sizeof(float) * 4));
- double needMem = needMemAll / double(comm->getNumberOfProcesses());
-
- if (myid == 0)
- {
- UBLOG(logINFO, "Number of blocks = " << numberOfBlocks);
- UBLOG(logINFO, "Number of nodes = " << numberOfNodes);
- int minInitLevel = grid->getCoarsestInitializedLevel();
- int maxInitLevel = grid->getFinestInitializedLevel();
- for (int level = minInitLevel; level <= maxInitLevel; level++)
+ SPtr<Communicator> comm = MPICommunicator::getInstance();
+ int myid = comm->getProcessID();
+
+ ConfigurationFile config;
+ config.load(configname);
+
+ vector<int> blocknx = config.getVector<int>("blocknx");
+ vector<double> boundingBox = config.getVector<double>("boundingBox");
+ peMinOffset = config.getVector<double>("peMinOffset");
+ peMaxOffset = config.getVector<double>("peMaxOffset");
+ int endTime = config.getValue<int>("endTime");
+ double outTime = config.getValue<double>("outTime");
+ double availMem = config.getValue<double>("availMem");
+ double uLB = config.getValue<double>("uLB");
+ pathOut = config.getValue<string>("pathOut");
+ pathGeo = config.getValue<string>("pathGeo");
+ string michel = config.getValue<string>("michel");
+ string plexiglas = config.getValue<string>("plexiglas");
+ double sphereTime = config.getValue<double>("sphereTime");
+
+ //parameters
+ //string pathOut = "d:/temp/thermoplast3";
+ //string pathGeo = "d:/Projects/ThermoPlast/Geometrie";
+ int numOfThreads = 1;
+ //int blocknx[3] ={ 10,10,10 };
+ //double endTime = 1000000;
+ //double outTime = 300;
+ //double availMem = 8e9;
+ double deltax = 1;
+ double rhoLB = 0.0;
+ //double uLB = 0.1;
+ double radius = 5;
+ double Re = 300;
+ double nuLB = (uLB*2.0*radius)/Re;
+
+
+
+ //geometry definition
+
+ //simulation bounding box
+ g_minX1 = boundingBox[0];
+ g_minX2 = boundingBox[1];
+ g_minX3 = boundingBox[2];
+
+ g_maxX1 = boundingBox[3];
+ g_maxX2 = boundingBox[4];
+ g_maxX3 = boundingBox[5];
+
+ double blockLength = blocknx[0]*deltax;
+
+ SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3));
+ if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), pathOut + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance());
+
+ //box
+ SPtr<GbObject3D> box(new GbCuboid3D(g_minX1-blockLength, g_minX2, g_minX3, g_maxX1+blockLength, g_maxX2, g_maxX3));
+ GbSystem3D::writeGeoObject(box.get(), pathOut + "/geo/box", WbWriterVtkXmlBinary::getInstance());
+
+ //michel
+ if (myid==0) UBLOG(logINFO, "Read michelGeo:start");
+ SPtr<GbTriFaceMesh3D> michelGeo = SPtr<GbTriFaceMesh3D>(GbTriFaceMesh3DCreator::getInstance()->readMeshFromSTLFile2(pathGeo+michel, "michelGeo", GbTriFaceMesh3D::KDTREE_SAHPLIT, false));
+ if (myid==0) UBLOG(logINFO, "Read michelGeo:end");
+ if (myid==0) GbSystem3D::writeGeoObject(michelGeo.get(), pathOut+"/geo/michelGeo", WbWriterVtkXmlBinary::getInstance());
+
+ //plexiglas
+ if (myid==0) UBLOG(logINFO, "Read plexiglasGeo:start");
+ SPtr<GbTriFaceMesh3D> plexiglasGeo = SPtr<GbTriFaceMesh3D>(GbTriFaceMesh3DCreator::getInstance()->readMeshFromSTLFile2(pathGeo+plexiglas, "plexiglasGeo", GbTriFaceMesh3D::KDTREE_SAHPLIT, false));
+ if (myid==0) UBLOG(logINFO, "Read plexiglasGeo:end");
+ if (myid==0) GbSystem3D::writeGeoObject(plexiglasGeo.get(), pathOut+"/geo/plexiglasGeo", WbWriterVtkXmlBinary::getInstance());
+
+ if (myid == 0)
{
- int nobl = grid->getNumberOfBlocks(level);
- UBLOG(logINFO, "Number of blocks for level " << level << " = " << nobl);
- UBLOG(logINFO, "Number of nodes for level " << level << " = " << nobl*numberOfNodesPerBlock);
+ UBLOG(logINFO, "Parameters:");
+ UBLOG(logINFO, "* uLB = " << uLB);
+ UBLOG(logINFO, "* rhoLB = " << rhoLB);
+ UBLOG(logINFO, "* nuLB = " << nuLB);
+ UBLOG(logINFO, "* deltaX = " << deltax);
+ UBLOG(logINFO, "* radius = " << radius);
+ UBLOG(logINFO, "* Re = " << Re);
+ UBLOG(logINFO, "* number of threads = "<<numOfThreads);
+ UBLOG(logINFO, "* number of processes = "<<comm->getNumberOfProcesses());
+ UBLOG(logINFO, "* path = "<<pathOut);
+ UBLOG(logINFO, "Preprocess - start");
}
- UBLOG(logINFO, "Necessary memory = " << needMemAll << " bytes");
- UBLOG(logINFO, "Necessary memory per process = " << needMem << " bytes");
- UBLOG(logINFO, "Available memory per process = " << availMem << " bytes");
- }
-
- //LBM kernel definition
- SPtr<LBMKernel> kernel;
- kernel = SPtr<LBMKernel>(new IncompressibleCumulantLBMKernel());
- SPtr<BCProcessor> bcProc(new BCProcessor());
- kernel->setBCProcessor(bcProc);
-
- //set forcing
- //mu::Parser fctForcingX1;
- //fctForcingX1.SetExpr("Fx1");
- //fctForcingX1.DefineConst("Fx1", 9e-7);
- //kernel->setWithForcing(true);
- //kernel->setForcingX1(fctForcingX1);
-
- //create LBM kernel
- SetKernelBlockVisitor kernelVisitor(kernel, nuLB, availMem, needMem);
- grid->accept(kernelVisitor);
-
- //set boundary conditions for nodes
- //michelInt->initInteractor();
- //plexiglasInt->initInteractor();
-
- intHelper.setBC();
- grid->accept(*bcVisitor.get());
-
- //initialization of distributions
- InitDistributionsBlockVisitor initVisitor;
- initVisitor.setVx1(uLB);
- grid->accept(initVisitor);
-
- //set connectors
- InterpolationProcessorPtr iProcessor(new IncompressibleOffsetInterpolationProcessor());
- SetConnectorsBlockVisitor setConnsVisitor(comm, true, D3Q27System::ENDDIR, nuLB, iProcessor);
- grid->accept(setConnsVisitor);
-
- //write data for visualization of boundary conditions
- {
- SPtr<UbScheduler> geoSch(new UbScheduler(1));
- WriteBoundaryConditionsCoProcessor ppgeo(grid, geoSch, pathOut, WbWriterVtkXmlBinary::getInstance(), comm);
- ppgeo.process(0);
- WriteMacroscopicQuantitiesCoProcessor ppInit(grid, geoSch, pathOut, WbWriterVtkXmlBinary::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm);
- ppInit.process(0);
- }
-
- if (myid == 0) UBLOG(logINFO, "Preprocess - end");
-
- //write data for visualization of macroscopic quantities
- SPtr<UbScheduler> visSch(new UbScheduler(outTime));
- SPtr<WriteMacroscopicQuantitiesCoProcessor> writeMQCoProcessor(new WriteMacroscopicQuantitiesCoProcessor(grid, visSch, pathOut,
- WbWriterVtkXmlBinary::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm));
-
- SPtr<WriteBoundaryConditionsCoProcessor> writeBCCoProcessor(new WriteBoundaryConditionsCoProcessor(grid, visSch, pathOut,
- WbWriterVtkXmlBinary::getInstance(), comm));
-
- SPtr<WriteDemObjectsCoProcessor> writeDemObjectsCoProcessor(new WriteDemObjectsCoProcessor(grid, visSch, pathOut, SPtr<WbWriter>(WbWriterVtkXmlBinary::getInstance()), demCoProcessor, comm));
- writeDemObjectsCoProcessor->process(0);
-
- //performance control
- SPtr<UbScheduler> nupsSch(new UbScheduler(10, 30, 100));
- SPtr<NUPSCounterCoProcessor> npr(new NUPSCounterCoProcessor(grid, nupsSch, numOfThreads, comm));
-
- //////////////////////////////////////////////////////////////////////////
- //generating spheres
- double d = 2.0*radius;
SPtr<UbScheduler> sphereScheduler(new UbScheduler(sphereTime));
Vector3D origin(g_minX1+60.0+4.0+radius, geoInflow1->getX2Minimum()+4.0*d, geoInflow1->getX3Minimum()+2.0*d);
- //std::array<double, 6> AABB={origin[0]-radius-1,origin[1]-radius,origin[2]-radius-1,origin[0]+radius+1, origin[1]+2.0*d+radius+1, origin[2]+2.0*d+radius+1};
- //SPtr<GbObject3D> boxAABB(new GbCuboid3D(AABB[0],AABB[1],AABB[2],AABB[3],AABB[4],AABB[5]));
- //GbSystem3D::writeGeoObject(boxAABB.get(), pathOut + "/geo/boxAABB", WbWriterVtkXmlBinary::getInstance());
- SPtr<CreateDemObjectsCoProcessor> createSphereCoProcessor(new CreateDemObjectsCoProcessor(grid,sphereScheduler,demCoProcessor,sphereMaterial,Vector3D(uLB, 0.0, 0.0)));
- //spheres
- for (int x3 = 0; x3 < 6; x3++)
for (int x2 = 0; x2 < 4; x2++)
for (int x1 = 0; x1 < 1; x1++)
+ //Grid definition
+ SPtr<Grid3D> grid(new Grid3D(comm));
+ grid->setDeltaX(deltax);
+ grid->setBlockNX(blocknx[0], blocknx[1], blocknx[2]);
+ grid->setPeriodicX1(false);
+ grid->setPeriodicX2(false);
+ grid->setPeriodicX3(false);
+
+ //blocks generating
+ GenBlocksGridVisitor genBlocks(gridCube);
+ grid->accept(genBlocks);
+
+ //boundary conditions definition
+ //boundary conditions adapters
+ //////////////////////////////////////////////////////////////////////////////
+ SPtr<BCAdapter> noSlipBCAdapter(new NoSlipBCAdapter());
+ noSlipBCAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NoSlipBCAlgorithm()));
+
+ mu::Parser fct;
+ fct.SetExpr("U");
+ fct.DefineConst("U", uLB);
+ SPtr<BCAdapter> inflowAdapter(new VelocityBCAdapter(true, false, false, fct, 0, BCFunction::INFCONST));
+ inflowAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityBCAlgorithm()));
+
+ SPtr<BCAdapter> outflowAdapter(new DensityBCAdapter(rhoLB));
+ outflowAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NonEqDensityBCAlgorithm()));
+ //outflowAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NonReflectingOutflowBCAlgorithm()));
+
+ //sphere
+ mu::Parser fct2;
+ fct2.SetExpr("U");
+ fct2.DefineConst("U", 0.0);
+ SPtr<BCAdapter> velocityBcParticleAdapter(new VelocityBCAdapter(true, false, false, fct2, 0, BCFunction::INFCONST));
+ velocityBcParticleAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityWithDensityBCAlgorithm()));
+
+ //inflow
+ GbCuboid3DPtr geoInflow1(new GbCuboid3D(g_minX1-blockLength, g_maxX2-120.0, g_minX3+190.0, g_minX1+1, g_maxX2+20.0, g_minX3+320.0));
+ if (myid == 0) GbSystem3D::writeGeoObject(geoInflow1.get(), pathOut + "/geo/geoInflow1", WbWriterVtkXmlASCII::getInstance());
+
+ GbCuboid3DPtr geoInflow2(new GbCuboid3D(g_minX1-blockLength, g_minX2 - blockLength, g_minX3 - blockLength, g_minX1, g_maxX2 + blockLength, g_maxX3 + blockLength));
+ if (myid == 0) GbSystem3D::writeGeoObject(geoInflow2.get(), pathOut + "/geo/geoInflow2", WbWriterVtkXmlASCII::getInstance());
+
+ GbCuboid3DPtr geoInflow3(new GbCuboid3D(g_minX1-blockLength, g_minX2-radius, g_maxX3-4.0*radius-1, g_minX1+1, g_maxX2+radius, g_maxX3-4.0*radius));
+ if (myid == 0) GbSystem3D::writeGeoObject(geoInflow3.get(), pathOut + "/geo/geoInflow3", WbWriterVtkXmlASCII::getInstance());
+
+ GbCuboid3DPtr geoInflow4(new GbCuboid3D(g_minX1-blockLength, g_minX2+4.0*radius, g_maxX3-4.0*radius-1.0, g_minX1+1, g_minX2+4.0*radius+1.0, g_maxX3+radius));
+ if (myid == 0) GbSystem3D::writeGeoObject(geoInflow4.get(), pathOut + "/geo/geoInflow4", WbWriterVtkXmlASCII::getInstance());
+
+ //outflow
+ GbCuboid3DPtr geoOutflow(new GbCuboid3D(g_maxX1, g_minX2 - blockLength, g_minX3 - blockLength, g_maxX1 + blockLength, g_maxX2 + blockLength, g_maxX3 + blockLength));
+ if (myid == 0) GbSystem3D::writeGeoObject(geoOutflow.get(), pathOut + "/geo/geoOutflow", WbWriterVtkXmlASCII::getInstance());
+
+ //boundary conditions visitor
+ SPtr<BoundaryConditionsBlockVisitor> bcVisitor(new BoundaryConditionsBlockVisitor());
+ bcVisitor->addBC(noSlipBCAdapter);
+ bcVisitor->addBC(inflowAdapter);
+ bcVisitor->addBC(outflowAdapter);
+ bcVisitor->addBC(velocityBcParticleAdapter);
+ //////////////////////////////////////////////////////////////////////////////////
+
+ //set boundary conditions for blocks and create process decomposition for MPI
+ SPtr<D3Q27Interactor> boxInt(new D3Q27Interactor(box, grid, noSlipBCAdapter, Interactor3D::INVERSESOLID));
+
+ //sphere bc object
+ //SPtr<MovableObjectInteractor> sphereInt1;
+ //const std::shared_ptr<Reconstructor> velocityReconstructor = std::make_shared<VelocityBcReconstructor>();
+ //const std::shared_ptr<Reconstructor> extrapolationReconstructor = std::make_shared<ExtrapolationReconstructor>(velocityReconstructor);
+ //sphereInt1 = SPtr<MovableObjectInteractor>(new MovableObjectInteractor(sphere1, grid, velocityBcParticleAdapter, Interactor3D::SOLID, extrapolationReconstructor, State::UNPIN));
+
+ //inflow
+ SPtr<D3Q27Interactor> inflowInt1 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow1, grid, inflowAdapter, Interactor3D::SOLID));
+ SPtr<D3Q27Interactor> inflowInt2 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow2, grid, outflowAdapter, Interactor3D::SOLID));
+
+ SPtr<D3Q27Interactor> inflowInt3 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow3, grid, noSlipBCAdapter, Interactor3D::SOLID));
+ SPtr<D3Q27Interactor> inflowInt4 = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow4, grid, noSlipBCAdapter, Interactor3D::SOLID));
+
+ //outflow
+ SPtr<D3Q27Interactor> outflowInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow, grid, outflowAdapter, Interactor3D::SOLID));
+
+ //michel
+ SPtr<Interactor3D> michelInt = SPtr<D3Q27TriFaceMeshInteractor>(new D3Q27TriFaceMeshInteractor(michelGeo, grid, noSlipBCAdapter, Interactor3D::SOLID));
+
+ //plexiglas
+ SPtr<Interactor3D> plexiglasInt = SPtr<D3Q27TriFaceMeshInteractor>(new D3Q27TriFaceMeshInteractor(plexiglasGeo, grid, noSlipBCAdapter, Interactor3D::SOLID));
+
+ //PE
+ double refLengthLb = radius*2.0;
+ double refLengthWorld = refLengthLb * deltax;
+ const std::shared_ptr<LBMUnitConverter> lbmUnitConverter = std::make_shared<LBMUnitConverter>(refLengthWorld, LBMUnitConverter::WORLD_MATERIAL::AIR_20C, refLengthLb);
+ if (myid == 0) std::cout << lbmUnitConverter->toString() << std::endl;
+ double rhoSphere = 915 * lbmUnitConverter->getFactorDensityWToLb(); // kg/m^3
+ if (myid == 0) UBLOG(logINFO, "rhoSphere = "<<rhoSphere);
+ SPtr<PhysicsEngineMaterialAdapter> sphereMaterial(new PePhysicsEngineMaterialAdapter("Polypropylen", rhoSphere, 0, 0.15, 0.1, 0.45, 0.5, 1, 0, 0));
+ const int timestep = 2;
+ const SPtr<UbScheduler> peScheduler(new UbScheduler(timestep));
+ int maxpeIterations = 10;//endTime/2;
+ SPtr<DemCoProcessor> demCoProcessor = makePeCoProcessor(grid, comm, peScheduler, lbmUnitConverter, maxpeIterations);
+ //demCoProcessor->addInteractor(sphereInt1, sphereMaterial, Vector3D(0.0, 0.0, 0.0));
+ //demCoProcessor->distributeIDs();
+ demCoProcessor->setBlockVisitor(bcVisitor);
+ //double g = 9.81 * lbmUnitConverter->getFactorAccWToLb();
+ //double rhoFluid = lbmUnitConverter->getFactorDensityWToLb() * 830.0; // 1 // kg/m^3
+ //////////////////////////////////////////////////////////////////////////
+ SPtr<Grid3DVisitor> peVisitor(new PePartitioningGridVisitor(comm, demCoProcessor));
+ //grid->accept(peVisitor);
+
+ //SetSolidBlocksBlockVisitor solidBlocksVisitor1(michelInt);
+ //grid->accept(solidBlocksVisitor1);
+ //SetSolidBlocksBlockVisitor solidBlocksVisitor2(plexiglasInt);
+ //grid->accept(solidBlocksVisitor2);
+
+ //SetBcBlocksBlockVisitor bcBlocksVisitor1(michelInt);
+ //grid->accept(bcBlocksVisitor1);
+ //SetBcBlocksBlockVisitor bcBlocksVisitor2(plexiglasInt);
+ //grid->accept(bcBlocksVisitor2);
+
+ InteractorsHelper intHelper(grid, peVisitor);
+ intHelper.addInteractor(boxInt);
+ //intHelper.addInteractor(sphereInt1);
+ //intHelper.addInteractor(inflowInt2);
+ //intHelper.addInteractor(michelInt);
+ //intHelper.addInteractor(plexiglasInt);
+ intHelper.addInteractor(inflowInt1);
+ intHelper.addInteractor(outflowInt);
+ intHelper.addInteractor(inflowInt2);
+ //intHelper.addInteractor(inflowInt3);
+ //intHelper.addInteractor(inflowInt4);
+
+ //intHelper.addInteractor(plexiglasInt);
+
+ intHelper.selectBlocks();
+
+ //write data for visualization of block grid
+ SPtr<CoProcessor> ppblocks(new WriteBlocksCoProcessor(grid, SPtr<UbScheduler>(new UbScheduler(1)), pathOut, WbWriterVtkXmlBinary::getInstance(), comm));
+ ppblocks->process(0);
+ ppblocks.reset();
+
+ unsigned long long numberOfBlocks = (unsigned long long)grid->getNumberOfBlocks();
+ int ghostLayer = 3;
+ unsigned long long numberOfNodesPerBlock = (unsigned long long)(blocknx[0])* (unsigned long long)(blocknx[1])* (unsigned long long)(blocknx[2]);
+ unsigned long long numberOfNodes = numberOfBlocks * numberOfNodesPerBlock;
+ unsigned long long numberOfNodesPerBlockWithGhostLayer = numberOfBlocks * (blocknx[0] + ghostLayer) * (blocknx[1] + ghostLayer) * (blocknx[2] + ghostLayer);
+ double needMemAll = double(numberOfNodesPerBlockWithGhostLayer*(27 * sizeof(double) + sizeof(int) + sizeof(float) * 4));
+ double needMem = needMemAll / double(comm->getNumberOfProcesses());
+
+ if (myid == 0)
+ {
+ UBLOG(logINFO, "Number of blocks = " << numberOfBlocks);
+ UBLOG(logINFO, "Number of nodes = " << numberOfNodes);
+ int minInitLevel = grid->getCoarsestInitializedLevel();
+ int maxInitLevel = grid->getFinestInitializedLevel();
+ for (int level = minInitLevel; level <= maxInitLevel; level++)
{
- //SPtr<GbObject3D> sphere(new GbSphere3D(origin[0]+x1*d, origin[1]+x2*2.0*d, origin[2]+x3*2.0*d, radius));
- SPtr<GbObject3D> sphere(new GbSphere3D(origin[0]+x1*d, origin[1]+x2*1.5*d, origin[2]+x3*1.5*d, radius));
- if (myid == 0) GbSystem3D::writeGeoObject(sphere.get(), pathOut + "/geo/sphere"+UbSystem::toString(x1)+UbSystem::toString(x2)+UbSystem::toString(x3), WbWriterVtkXmlASCII::getInstance());
- createSphereCoProcessor->addGeoObject(sphere);
+ int nobl = grid->getNumberOfBlocks(level);
+ UBLOG(logINFO, "Number of blocks for level " << level << " = " << nobl);
+ UBLOG(logINFO, "Number of nodes for level " << level << " = " << nobl*numberOfNodesPerBlock);
}
+ UBLOG(logINFO, "Necessary memory = " << needMemAll << " bytes");
+ UBLOG(logINFO, "Necessary memory per process = " << needMem << " bytes");
+ UBLOG(logINFO, "Available memory per process = " << availMem << " bytes");
+ }
- createSphereCoProcessor->process(0);
-
- //start simulation
- omp_set_num_threads(numOfThreads);
- SPtr<UbScheduler> stepGhostLayer(visSch);
- SPtr<Calculator> calculator(new BasicCalculator(grid, peScheduler, endTime));
- calculator->addCoProcessor(npr);
-
- calculator->addCoProcessor(createSphereCoProcessor);
- calculator->addCoProcessor(demCoProcessor);
-
- calculator->addCoProcessor(writeBCCoProcessor);
- calculator->addCoProcessor(writeDemObjectsCoProcessor);
- calculator->addCoProcessor(writeMQCoProcessor);
-
- if (myid == 0) UBLOG(logINFO, "Simulation-start");
- calculator->calculate();
- if (myid == 0) UBLOG(logINFO, "Simulation-end");
+ //LBM kernel definition
+ SPtr<LBMKernel> kernel;
+ kernel = SPtr<LBMKernel>(new IncompressibleCumulantLBMKernel());
+ SPtr<BCProcessor> bcProc(new BCProcessor());
+ kernel->setBCProcessor(bcProc);
+
+ //set forcing
+ //mu::Parser fctForcingX1;
+ //fctForcingX1.SetExpr("Fx1");
+ //fctForcingX1.DefineConst("Fx1", 9e-7);
+ //kernel->setWithForcing(true);
+ //kernel->setForcingX1(fctForcingX1);
+
+ //create LBM kernel
+ SetKernelBlockVisitor kernelVisitor(kernel, nuLB, availMem, needMem);
+ grid->accept(kernelVisitor);
+
+ //set boundary conditions for nodes
+ //michelInt->initInteractor();
+ //plexiglasInt->initInteractor();
+
+ intHelper.setBC();
+ grid->accept(*bcVisitor.get());
+
+ //initialization of distributions
+ InitDistributionsBlockVisitor initVisitor;
+ //initVisitor.setVx1(uLB);
+ grid->accept(initVisitor);
+
+ //set connectors
+ InterpolationProcessorPtr iProcessor(new IncompressibleOffsetInterpolationProcessor());
+ SetConnectorsBlockVisitor setConnsVisitor(comm, true, D3Q27System::ENDDIR, nuLB, iProcessor);
+ grid->accept(setConnsVisitor);
+
+ //write data for visualization of boundary conditions
+ {
+ SPtr<UbScheduler> geoSch(new UbScheduler(1));
+ WriteBoundaryConditionsCoProcessor ppgeo(grid, geoSch, pathOut, WbWriterVtkXmlBinary::getInstance(), comm);
+ ppgeo.process(0);
+ WriteMacroscopicQuantitiesCoProcessor ppInit(grid, geoSch, pathOut, WbWriterVtkXmlBinary::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm);
+ ppInit.process(0);
}
- catch (std::exception& e)
{
UBLOG(logERROR, e.what());
}
catch (std::string& s)
{
UBLOG(logERROR, s);
}
catch (...)
{
UBLOG(logERROR, "unknown exception");
}
+
+ if (myid == 0) UBLOG(logINFO, "Preprocess - end");
+
+ //write data for visualization of macroscopic quantities
+ SPtr<UbScheduler> visSch(new UbScheduler(outTime));
+ SPtr<WriteMacroscopicQuantitiesCoProcessor> writeMQCoProcessor(new WriteMacroscopicQuantitiesCoProcessor(grid, visSch, pathOut,
+ WbWriterVtkXmlBinary::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm));
+
+ SPtr<WriteBoundaryConditionsCoProcessor> writeBCCoProcessor(new WriteBoundaryConditionsCoProcessor(grid, visSch, pathOut,
+ WbWriterVtkXmlBinary::getInstance(), comm));
+
+ SPtr<WriteDemObjectsCoProcessor> writeDemObjectsCoProcessor(new WriteDemObjectsCoProcessor(grid, visSch, pathOut, SPtr<WbWriter>(WbWriterVtkXmlBinary::getInstance()), demCoProcessor, comm));
+ writeDemObjectsCoProcessor->process(0);
+
+ //performance control
+ SPtr<UbScheduler> nupsSch(new UbScheduler(10, 30, 100));
+ SPtr<NUPSCounterCoProcessor> npr(new NUPSCounterCoProcessor(grid, nupsSch, numOfThreads, comm));
+
+ //////////////////////////////////////////////////////////////////////////
+ //generating spheres
+ double d = 2.0*radius;
+ SPtr<UbScheduler> sphereScheduler(new UbScheduler(sphereTime));
+ Vector3D origin(g_minX1+60.0+4.0+radius, geoInflow1->getX2Minimum()+4.0*d, geoInflow1->getX3Minimum()+2.0*d);
+ //std::array<double, 6> AABB={origin[0]-radius-1,origin[1]-radius,origin[2]-radius-1,origin[0]+radius+1, origin[1]+2.0*d+radius+1, origin[2]+2.0*d+radius+1};
+ //SPtr<GbObject3D> boxAABB(new GbCuboid3D(AABB[0],AABB[1],AABB[2],AABB[3],AABB[4],AABB[5]));
+ //GbSystem3D::writeGeoObject(boxAABB.get(), pathOut + "/geo/boxAABB", WbWriterVtkXmlBinary::getInstance());
+ SPtr<CreateDemObjectsCoProcessor> createSphereCoProcessor(new CreateDemObjectsCoProcessor(grid, sphereScheduler, demCoProcessor, sphereMaterial));
+ //spheres
+ //for (int x3 = 0; x3 < 6; x3++)
+ // for (int x2 = 0; x2 < 4; x2++)
+ // for (int x1 = 0; x1 < 1; x1++)
+ // {
+ // //SPtr<GbObject3D> sphere(new GbSphere3D(origin[0]+x1*d, origin[1]+x2*2.0*d, origin[2]+x3*2.0*d, radius));
+ // SPtr<GbObject3D> sphere(new GbSphere3D(origin[0]+x1*d, origin[1]+x2*1.5*d, origin[2]+x3*1.5*d, radius));
+ // if (myid == 0) GbSystem3D::writeGeoObject(sphere.get(), pathOut + "/geo/sphere"+UbSystem::toString(x1)+UbSystem::toString(x2)+UbSystem::toString(x3), WbWriterVtkXmlASCII::getInstance());
+ // createSphereCoProcessor->addGeoObject(sphere, Vector3D(uLB, 0.0, 0.0));
+ // }
+
+ //createSphereCoProcessor->process(0);
+
+ SPtr<MPIIORestartCoProcessor> restartCoProcessor(new MPIIORestartCoProcessor(grid, sphereScheduler, pathOut, comm));
+ restartCoProcessor->setLBMKernel(kernel);
+ restartCoProcessor->setBCProcessor(bcProc);
+
+ SPtr<RestartDemObjectsCoProcessor> restartDemObjectsCoProcessor(new RestartDemObjectsCoProcessor(grid, sphereScheduler, pathOut, demCoProcessor, createSphereCoProcessor, radius, comm));
+
+ restartCoProcessor->restart(10);
+ restartDemObjectsCoProcessor->read(10);
+
+
+
+ // for (int x3 = 0; x3 < 6; x3++)
+ //for (int x2 = 0; x2 < 4; x2++)
+ // for (int x1 = 0; x1 < 1; x1++)
+ // {
+ // //SPtr<GbObject3D> sphere(new GbSphere3D(origin[0]+x1*d, origin[1]+x2*2.0*d, origin[2]+x3*2.0*d, radius));
+ // SPtr<GbObject3D> sphere(new GbSphere3D(origin[0]+x1*d, origin[1]+x2*1.5*d, origin[2]+x3*1.5*d, radius));
+ // if (myid == 0) GbSystem3D::writeGeoObject(sphere.get(), pathOut + "/geo/sphere"+UbSystem::toString(x1)+UbSystem::toString(x2)+UbSystem::toString(x3), WbWriterVtkXmlASCII::getInstance());
+ // createSphereCoProcessor->addGeoObject(sphere, Vector3D(uLB, 0.0, 0.0));
+ // }
+
+ //start simulation
+ omp_set_num_threads(numOfThreads);
+ SPtr<UbScheduler> stepGhostLayer(visSch);
+ SPtr<Calculator> calculator(new BasicCalculator(grid, peScheduler, endTime));
+ calculator->addCoProcessor(npr);
+
+ calculator->addCoProcessor(createSphereCoProcessor);
+ calculator->addCoProcessor(demCoProcessor);
+ //calculator->addCoProcessor(writeBCCoProcessor);
+ calculator->addCoProcessor(writeDemObjectsCoProcessor);
+ calculator->addCoProcessor(writeMQCoProcessor);
+ calculator->addCoProcessor(restartDemObjectsCoProcessor);
+ calculator->addCoProcessor(restartCoProcessor);
+
+ if (myid == 0) UBLOG(logINFO, "Simulation-start");
+ calculator->calculate();
+ if (myid == 0) UBLOG(logINFO, "Simulation-end");
+
+ }
+ catch (std::exception& e)
+ {
+ UBLOG(logERROR, e.what());
+ }
+ catch (std::string& s)
+ {
+ UBLOG(logERROR, s);
+ }
+ catch (...)
+ {
+ UBLOG(logERROR, "unknown exception");
+ }
}
@@ -446,20 +485,31 @@ int main(int argc, char* argv[])
//try
//{
//Sleep(30000);
- walberla::Environment env(argc, argv);
+ walberla::Environment env(argc, argv);
- if (argv!=NULL)
+ if (argv!=NULL)
+ {
+ if (argv[1]!=NULL)
{
- if (argv[1]!=NULL)
- {
- thermoplast(string(argv[1]));
- }
- else
- {
- cout<<"Configuration file must be set!: "<<argv[0]<<" <config file>"<<endl<<std::flush;
- }
+ thermoplast(string(argv[1]));
}
- return 0;
+ else
+ {
+ cout<<"Configuration file must be set!: "<<argv[0]<<" <config file>"<<endl<<std::flush;
+ }
+ }
+ return 0;
+ //}
+ //catch (std::exception& e)
+ //{
+ // UBLOG(logERROR, e.what());
+ //}
+ //catch (std::string& s)
+ //{
+ // UBLOG(logERROR, s);
+ //}
+ //catch (...)
+ //{
+ // UBLOG(logERROR, "unknown exception");
//}
- //catch (std::exception& e)
//{
// UBLOG(logERROR, e.what());
//}
//catch (std::string& s)
//{
// UBLOG(logERROR, s);
//}
//catch (...)
//{
// UBLOG(logERROR, "unknown exception");
//}
}
diff --git a/source/DemCoupling/CreateDemObjectsCoProcessor.cpp b/source/DemCoupling/CreateDemObjectsCoProcessor.cpp
index 28311218eeaa85cb8440d008c8e88ec2adf3bbb8..cab9e39d717864f63482b0d2c4ceb33e91c70c69 100644
--- a/source/DemCoupling/CreateDemObjectsCoProcessor.cpp
+++ b/source/DemCoupling/CreateDemObjectsCoProcessor.cpp
@@ -17,52 +17,90 @@
#include "SetBcBlocksBlockVisitor.h"
#include "Grid3D.h"
-CreateDemObjectsCoProcessor::CreateDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, SPtr<DemCoProcessor> demCoProcessor, SPtr<PhysicsEngineMaterialAdapter> demObjectMaterial, Vector3D initalVelocity) :
+CreateDemObjectsCoProcessor::CreateDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, SPtr<DemCoProcessor> demCoProcessor, SPtr<PhysicsEngineMaterialAdapter> demObjectMaterial) :
CoProcessor(grid, s),
demCoProcessor(demCoProcessor),
demObjectMaterial(demObjectMaterial),
initalVelocity(initalVelocity)
{
+ mu::Parser fct;
+ fct.SetExpr("U");
+ fct.DefineConst("U", 0.0);
+ velocityBcParticleAdapter = SPtr<BCAdapter>(new VelocityBCAdapter(true, false, false, fct, 0, BCFunction::INFCONST));
+ velocityBcParticleAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityWithDensityBCAlgorithm()));
+ //const std::shared_ptr<Reconstructor> velocityReconstructor(new VelocityBcReconstructor());
+ std::shared_ptr<Reconstructor> equilibriumReconstructor(new EquilibriumReconstructor());
+ //const std::shared_ptr<Reconstructor> lbmReconstructor(new LBMReconstructor(false));
+ extrapolationReconstructor = SPtr<Reconstructor>(new ExtrapolationReconstructor(equilibriumReconstructor));
}
//////////////////////////////////////////////////////////////////////////
void CreateDemObjectsCoProcessor::process(double step)
{
if (scheduler->isDue(step))
{
- mu::Parser fct2;
- fct2.SetExpr("U");
- fct2.DefineConst("U", 0.0);
- SPtr<BCAdapter> velocityBcParticleAdapter(new VelocityBCAdapter(true, false, false, fct2, 0, BCFunction::INFCONST));
- velocityBcParticleAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityWithDensityBCAlgorithm()));
+ createGeoObjects();
+
+ //mu::Parser fct2;
+ //fct2.SetExpr("U");
+ //fct2.DefineConst("U", 0.0);
+ //SPtr<BCAdapter> velocityBcParticleAdapter(new VelocityBCAdapter(true, false, false, fct2, 0, BCFunction::INFCONST));
+ //velocityBcParticleAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new VelocityWithDensityBCAlgorithm()));
//velocityBcParticleAdapter->setBcAlgorithm(SPtr<BCAlgorithm>(new NoSlipBCAlgorithm()));
- //const std::shared_ptr<Reconstructor> velocityReconstructor(new VelocityBcReconstructor());
- const std::shared_ptr<Reconstructor> equilibriumReconstructor(new EquilibriumReconstructor());
- //const std::shared_ptr<Reconstructor> lbmReconstructor(new LBMReconstructor(false));
- const std::shared_ptr<Reconstructor> extrapolationReconstructor(new ExtrapolationReconstructor(equilibriumReconstructor));
+ ////const std::shared_ptr<Reconstructor> velocityReconstructor(new VelocityBcReconstructor());
+ //const std::shared_ptr<Reconstructor> equilibriumReconstructor(new EquilibriumReconstructor());
+ ////const std::shared_ptr<Reconstructor> lbmReconstructor(new LBMReconstructor(false));
+ //const std::shared_ptr<Reconstructor> extrapolationReconstructor(new ExtrapolationReconstructor(equilibriumReconstructor));
- for(SPtr<GbObject3D> proto : geoObjectPrototypeVector)
- {
- std::array<double, 6> AABB = {proto->getX1Minimum(),proto->getX2Minimum(),proto->getX3Minimum(),proto->getX1Maximum(),proto->getX2Maximum(),proto->getX3Maximum()};
- //UBLOG(logINFO, "demCoProcessor->isGeoObjectInAABB(AABB) = " << demCoProcessor->isGeoObjectInAABB(AABB));
- if (demCoProcessor->isDemObjectInAABB(AABB))
- {
- continue;
- }
- SPtr<GbObject3D> geoObject((GbObject3D*)(proto->clone()));
- SPtr<MovableObjectInteractor> geoObjectInt = SPtr<MovableObjectInteractor>(new MovableObjectInteractor(geoObject, grid, velocityBcParticleAdapter, Interactor3D::SOLID, extrapolationReconstructor, State::UNPIN));
- SetBcBlocksBlockVisitor setBcVisitor(geoObjectInt);
- grid->accept(setBcVisitor);
- geoObjectInt->initInteractor();
- demCoProcessor->addInteractor(geoObjectInt, demObjectMaterial, initalVelocity);
- //demCoProcessor->distributeIDs();
- }
+ //for(SPtr<GbObject3D> proto : geoObjectPrototypeVector)
+ //{
+ // std::array<double, 6> AABB = {proto->getX1Minimum(),proto->getX2Minimum(),proto->getX3Minimum(),proto->getX1Maximum(),proto->getX2Maximum(),proto->getX3Maximum()};
+ // //UBLOG(logINFO, "demCoProcessor->isGeoObjectInAABB(AABB) = " << demCoProcessor->isGeoObjectInAABB(AABB));
+ // if (demCoProcessor->isDemObjectInAABB(AABB))
+ // {
+ // continue;
+ // }
+ // SPtr<GbObject3D> geoObject((GbObject3D*)(proto->clone()));
+ // SPtr<MovableObjectInteractor> geoObjectInt = SPtr<MovableObjectInteractor>(new MovableObjectInteractor(geoObject, grid, velocityBcParticleAdapter, Interactor3D::SOLID, extrapolationReconstructor, State::UNPIN));
+ // SetBcBlocksBlockVisitor setBcVisitor(geoObjectInt);
+ // grid->accept(setBcVisitor);
+ // geoObjectInt->initInteractor();
+ // demCoProcessor->addInteractor(geoObjectInt, demObjectMaterial, initalVelocity);
+ //}
}
}
//////////////////////////////////////////////////////////////////////////
-void CreateDemObjectsCoProcessor::addGeoObject(SPtr<GbObject3D> geoObjectPrototype)
+void CreateDemObjectsCoProcessor::addGeoObject(SPtr<GbObject3D> geoObjectPrototype, Vector3D initalVelocity)
{
geoObjectPrototypeVector.push_back(geoObjectPrototype);
+ this->initalVelocity.push_back(initalVelocity);
+}
+
+void CreateDemObjectsCoProcessor::clearGeoObjects()
+{
+ geoObjectPrototypeVector.clear();
+ initalVelocity.clear();
+}
+
+void CreateDemObjectsCoProcessor::createGeoObjects()
+{
+ int size = geoObjectPrototypeVector.size();
+
+ for (int i = 0; i < size; i++)
+ {
+ std::array<double, 6> AABB ={ geoObjectPrototypeVector[i]->getX1Minimum(),geoObjectPrototypeVector[i]->getX2Minimum(),geoObjectPrototypeVector[i]->getX3Minimum(),geoObjectPrototypeVector[i]->getX1Maximum(),geoObjectPrototypeVector[i]->getX2Maximum(),geoObjectPrototypeVector[i]->getX3Maximum() };
+ //UBLOG(logINFO, "demCoProcessor->isGeoObjectInAABB(AABB) = " << demCoProcessor->isGeoObjectInAABB(AABB));
+ if (demCoProcessor->isDemObjectInAABB(AABB))
+ {
+ continue;
+ }
+ SPtr<GbObject3D> geoObject((GbObject3D*)(geoObjectPrototypeVector[i]->clone()));
+ SPtr<MovableObjectInteractor> geoObjectInt = SPtr<MovableObjectInteractor>(new MovableObjectInteractor(geoObject, grid, velocityBcParticleAdapter, Interactor3D::SOLID, extrapolationReconstructor, State::UNPIN));
+ SetBcBlocksBlockVisitor setBcVisitor(geoObjectInt);
+ grid->accept(setBcVisitor);
+ geoObjectInt->initInteractor();
+ demCoProcessor->addInteractor(geoObjectInt, demObjectMaterial, initalVelocity[i]);
+ }
}
diff --git a/source/DemCoupling/CreateDemObjectsCoProcessor.h b/source/DemCoupling/CreateDemObjectsCoProcessor.h
index 6b6a6ec15a6b7f3fd9b66660831f7a019a9508f3..01d5037dce165fbc9297d6e213867fdb7c66c4b4 100644
--- a/source/DemCoupling/CreateDemObjectsCoProcessor.h
+++ b/source/DemCoupling/CreateDemObjectsCoProcessor.h
@@ -10,20 +10,26 @@ class Grid3D;
class UbScheduler;
class DemCoProcessor;
class GbObject3D;
+class BCAdapter;
+class Reconstructor;
class PhysicsEngineMaterialAdapter;
class CreateDemObjectsCoProcessor : public CoProcessor
{
public:
- CreateDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, SPtr<DemCoProcessor> demCoProcessor, SPtr<PhysicsEngineMaterialAdapter> geoObjectMaterial, Vector3D initalVelocity);
+ CreateDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, SPtr<DemCoProcessor> demCoProcessor, SPtr<PhysicsEngineMaterialAdapter> geoObjectMaterial);
void process(double step) override;
- void addGeoObject( SPtr<GbObject3D> geoObjectPrototype);
+ void addGeoObject(SPtr<GbObject3D> geoObjectPrototype, Vector3D initalVelocity);
+ void clearGeoObjects();
+ void createGeoObjects();
protected:
private:
SPtr<DemCoProcessor> demCoProcessor;
std::vector< SPtr<GbObject3D> > geoObjectPrototypeVector;
SPtr<PhysicsEngineMaterialAdapter> demObjectMaterial;
- Vector3D initalVelocity;
+ std::vector<Vector3D> initalVelocity;
+ SPtr<BCAdapter> velocityBcParticleAdapter;
+ SPtr<Reconstructor> extrapolationReconstructor;
};
#endif // CreateSphereCoProcessor_h__
diff --git a/source/DemCoupling/DemCoProcessor.cpp b/source/DemCoupling/DemCoProcessor.cpp
index dedec708fcc829ea06c09968df42f98e89654739..ce4de51cdf582ee2fc8bcfefd2739f80da6f4f55 100644
--- a/source/DemCoupling/DemCoProcessor.cpp
+++ b/source/DemCoupling/DemCoProcessor.cpp
@@ -51,7 +51,7 @@ void DemCoProcessor::addInteractor(std::shared_ptr<MovableObjectInteractor> inte
{
physicsEngineGeometries.push_back(peGeometry);
}
- distributeIDs();
+ //distributeIDs();
}
@@ -266,6 +266,23 @@ void DemCoProcessor::addSurfaceTriangleSet(std::vector<UbTupleFloat3>& nodes, st
}
}
+void DemCoProcessor::getObjectsPropertiesVector(std::vector<double>& p)
+{
+ for (int i = 0; i < interactors.size(); i++)
+ {
+ if (std::dynamic_pointer_cast<PePhysicsEngineGeometryAdapter>(physicsEngineGeometries[i])->isActive())
+ {
+ p.push_back(interactors[i]->getGbObject3D()->getX1Centroid());
+ p.push_back(interactors[i]->getGbObject3D()->getX2Centroid());
+ p.push_back(interactors[i]->getGbObject3D()->getX3Centroid());
+ Vector3D v = physicsEngineGeometries[i]->getLinearVelocity();
+ p.push_back(v[0]);
+ p.push_back(v[1]);
+ p.push_back(v[2]);
+ }
+ }
+}
+
void DemCoProcessor::distributeIDs()
{
std::vector<int> peIDsSend;
@@ -295,7 +312,15 @@ void DemCoProcessor::distributeIDs()
for (int i = 0; i < interactors.size(); i++)
{
- std::dynamic_pointer_cast<PePhysicsEngineGeometryAdapter>(physicsEngineGeometries[i])->setId(idMap.find(interactors[i]->getID())->second);
+ std::map<int, int>::const_iterator it;
+ if ((it=idMap.find(interactors[i]->getID())) == idMap.end())
+ {
+ throw UbException(UB_EXARGS, "not valid ID!");
+ }
+
+
+ std::dynamic_pointer_cast<PePhysicsEngineGeometryAdapter>(physicsEngineGeometries[i])->setId(it->second);
+ //std::dynamic_pointer_cast<PePhysicsEngineGeometryAdapter>(physicsEngineGeometries[i])->setId(idMap.find(interactors[i]->getID())->second);
}
for (int i = 0; i < physicsEngineGeometries.size(); i++)
diff --git a/source/DemCoupling/DemCoProcessor.h b/source/DemCoupling/DemCoProcessor.h
index 8df8df76a7a22fdfa3dd50a03ddc05b51df0423e..aa9836113ceddc0227233118abf0c343fbb49282 100644
--- a/source/DemCoupling/DemCoProcessor.h
+++ b/source/DemCoupling/DemCoProcessor.h
@@ -40,6 +40,7 @@ public:
void setBlockVisitor(std::shared_ptr<BoundaryConditionsBlockVisitor> blockVisitor);
bool isDemObjectInAABB(std::array<double,6> AABB);
void addSurfaceTriangleSet(std::vector<UbTupleFloat3>& nodes, std::vector<UbTupleInt3>& triangles);
+ void getObjectsPropertiesVector(std::vector<double>& p);
private:
std::shared_ptr<PhysicsEngineGeometryAdapter> createPhysicsEngineGeometryAdapter(std::shared_ptr<MovableObjectInteractor> interactor, std::shared_ptr<PhysicsEngineMaterialAdapter> physicsEngineMaterial) const;
diff --git a/source/DemCoupling/DemCoupling.cmake b/source/DemCoupling/DemCoupling.cmake
index 12762704e3cdbeaf03d5feea5d260918cc2c4869..0fa6b211ba4b2281242eeb0f59843d1fa69fa9fd 100644
--- a/source/DemCoupling/DemCoupling.cmake
+++ b/source/DemCoupling/DemCoupling.cmake
@@ -21,3 +21,8 @@ SET(CAB_ADDITIONAL_LINK_LIBRARIES ${CAB_ADDITIONAL_LINK_LIBRARIES} ${LINK_LIBRAR
IF(${USE_GCC})
SET(CAB_ADDITIONAL_LINK_LIBRARIES ${CAB_ADDITIONAL_LINK_LIBRARIES} "stdc++fs")
ENDIF()
+
+IF(${USE_METIS})
+ SET(LINK_LIBRARY optimized ${METIS_RELEASE_LIBRARY} debug ${METIS_DEBUG_LIBRARY})
+ SET(CAB_ADDITIONAL_LINK_LIBRARIES ${CAB_ADDITIONAL_LINK_LIBRARIES} ${LINK_LIBRARY})
+ENDIF()
\ No newline at end of file
diff --git a/source/DemCoupling/RestartDemCouplingCoProcessor.cpp b/source/DemCoupling/RestartDemCouplingCoProcessor.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..9b657c21264922359d38d63b930b4aa6d354776c
--- /dev/null
+++ b/source/DemCoupling/RestartDemCouplingCoProcessor.cpp
@@ -0,0 +1,90 @@
+#include "RestartDemCouplingCoProcessor.h"
+
+#include "Vector3D.h"
+#include "Communicator.h"
+#include "UbScheduler.h"
+#include "Grid3D.h"
+#include "UbSystem.h"
+#include "GbSphere3D.h"
+#include "DemCoProcessor.h"
+#include "UbFileInputBinary.h"
+#include "UbFileOutputBinary.h"
+#include "CreateDemObjectsCoProcessor.h"
+
+RestartDemObjectsCoProcessor::RestartDemObjectsCoProcessor()
+{
+}
+
+RestartDemObjectsCoProcessor::RestartDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, const std::string & path, SPtr<DemCoProcessor> demCoProcessor, SPtr<CreateDemObjectsCoProcessor> createDemObjectsCoProcessor, double radius, SPtr<Communicator> comm) : CoProcessor(grid, s), path(path), demCoProcessor(demCoProcessor), createDemObjectsCoProcessor(createDemObjectsCoProcessor), radius(radius), comm(comm)
+{
+}
+
+void RestartDemObjectsCoProcessor::process(double step)
+{
+ if (scheduler->isDue(step))
+ {
+ int istep = static_cast<int>(step);
+
+ write(istep);
+
+ if (comm->isRoot())
+ UBLOG(logINFO, "RestartDemObjectsCoProcessor write step: " << istep);
+ }
+}
+
+void RestartDemObjectsCoProcessor::write(int step)
+{
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor write:start ");
+ std::vector<double> p;
+
+ demCoProcessor->getObjectsPropertiesVector(p);
+
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor size p: " << p.size());
+
+ std::vector<double> values = comm->gather(p);
+
+ if (comm->isRoot())
+ {
+ std::string subfolder = "dem_cp_"+UbSystem::toString(step);
+ std::string filePath = path+"/dem_cp/"+subfolder+"/dem_cp.bin";
+ UbFileOutputBinary fo(filePath);
+ fo.writeInteger((int)values.size());
+ fo.writeVector<double>(values);
+ UBLOG(logINFO, "RestartDemObjectsCoProcessor size: " << values.size());
+ }
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor write:stop ");
+}
+
+void RestartDemObjectsCoProcessor::read(int step)
+{
+ std::vector<double> p;
+
+ if (comm->isRoot())
+ {
+ std::string subfolder = "dem_cp_"+UbSystem::toString(step);
+ std::string filePath = path+"/dem_cp/"+subfolder+"/dem_cp.bin";
+ UbFileInputBinary fi(filePath);
+ int size = fi.readInteger();
+ p.resize(size);
+ fi.readVector<double>(p);
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor read size p: " << p.size());
+ }
+ comm->broadcast(p);
+
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor read size p broadcast: " << p.size());
+
+ createDemObjectsCoProcessor->clearGeoObjects();
+
+ int size = (int)p.size();
+
+ for (int i = 0; i < size; i += 6)
+ {
+ SPtr<GbObject3D> sphere(new GbSphere3D(p[i], p[i+1], p[i+2], radius));
+ createDemObjectsCoProcessor->addGeoObject(sphere, Vector3D(p[i+3], p[i+4], p[i+5]));
+ }
+
+ createDemObjectsCoProcessor->createGeoObjects();
+
+ createDemObjectsCoProcessor->clearGeoObjects();
+
+}
diff --git a/source/DemCoupling/RestartDemCouplingCoProcessor.h b/source/DemCoupling/RestartDemCouplingCoProcessor.h
new file mode 100644
index 0000000000000000000000000000000000000000..c44b2d7339653b7c85227fdb1290871297ca2864
--- /dev/null
+++ b/source/DemCoupling/RestartDemCouplingCoProcessor.h
@@ -0,0 +1,37 @@
+/*
+* Author: K. Kutscher
+* mail: kutscher@irmb.tu-bs.de
+*/
+#ifndef RestartDemObjectsCoProcessor_H
+#define RestartDemObjectsCoProcessor_H
+
+#include <PointerDefinitions.h>
+#include <string>
+#include <vector>
+
+#include "CoProcessor.h"
+
+class Communicator;
+class Grid3D;
+class UbScheduler;
+class DemCoProcessor;
+class CreateDemObjectsCoProcessor;
+
+class RestartDemObjectsCoProcessor : public CoProcessor
+{
+public:
+ RestartDemObjectsCoProcessor();
+ RestartDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, const std::string& path, SPtr<DemCoProcessor> demCoProcessor, SPtr<CreateDemObjectsCoProcessor> createDemObjectsCoProcessor, double radius, SPtr<Communicator> comm);
+ ~RestartDemObjectsCoProcessor() {}
+ void process(double step) override;
+ void write(int step);
+ void read(int step);
+
+private:
+ std::string path;
+ double radius;
+ SPtr<Communicator> comm;
+ SPtr<DemCoProcessor> demCoProcessor;
+ SPtr<CreateDemObjectsCoProcessor> createDemObjectsCoProcessor;
+};
+#endif
diff --git a/source/DemCoupling/RestartDemObjectsCoProcessor.cpp b/source/DemCoupling/RestartDemObjectsCoProcessor.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/source/DemCoupling/WriteDemObjectsCoProcessor.cpp b/source/DemCoupling/WriteDemObjectsCoProcessor.cpp
index 89672213b31b0434c194364ff55923ffa68bd2bd..682433a2e1dafd3bb40234b7ad9c23112d359083 100644
--- a/source/DemCoupling/WriteDemObjectsCoProcessor.cpp
+++ b/source/DemCoupling/WriteDemObjectsCoProcessor.cpp
@@ -1,6 +1,4 @@
#include "WriteDemObjectsCoProcessor.h"
-#include "LBMKernel.h"
-#include "BCProcessor.h"
#include "basics/writer/WbWriterVtkXmlBinary.h"
#include "basics/writer/WbWriterVtkXmlASCII.h"
@@ -23,7 +21,7 @@ WriteDemObjectsCoProcessor::WriteDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<U
demCoProcessor(demCoProcessor),
comm(comm)
{
- //this->path += "/geo/dem/objects";
+
}
//////////////////////////////////////////////////////////////////////////
void WriteDemObjectsCoProcessor::process(double step)
@@ -71,8 +69,5 @@ void WriteDemObjectsCoProcessor::process(double step)
}
UBLOG(logINFO, "WriteDemObjectsCoProcessor step: " << istep);
}
-
}
-
-
}
diff --git a/source/DemCoupling/physicsEngineAdapter/pe/PePhysicsEngineSolverAdapter.cpp b/source/DemCoupling/physicsEngineAdapter/pe/PePhysicsEngineSolverAdapter.cpp
index 0f15b8bae271da0944872ab7d62a7b47c2719f70..c85219e05d63a24920c00631a3b0ecd37af2c1df 100644
--- a/source/DemCoupling/physicsEngineAdapter/pe/PePhysicsEngineSolverAdapter.cpp
+++ b/source/DemCoupling/physicsEngineAdapter/pe/PePhysicsEngineSolverAdapter.cpp
@@ -10,6 +10,8 @@
#include "UbLogger.h"
#include <boost/tuple/tuple.hpp>
+#include <memory>
+
typedef boost::tuple<walberla::pe::Sphere, walberla::pe::Plane> BodyTypeTuple;
@@ -176,3 +178,30 @@ std::shared_ptr< walberla::blockforest::BlockForest > PePhysicsEngineSolverAdapt
{
return forest;
}
+
+void PePhysicsEngineSolverAdapter::saveToFile(const std::string & path)
+{
+ //forest->saveToFile(path+"SerializeDeserialize.sbf");
+ //forest->saveBlockData("SerializeDeserialize.dump", *storageId.get());
+}
+
+void PePhysicsEngineSolverAdapter::loadFromFile(const std::string & path)
+{
+ //forest = std::make_shared< walberla::blockforest::BlockForest >( walberla::uint_c( walberla::MPIManager::instance()->rank() ), path+"SerializeDeserialize.sbf", true, false );
+ //storageId = std::make_shared< walberla::domain_decomposition::BlockDataID >(forest->loadBlockData(path+"SerializeDeserialize.dump", walberla::pe::createStorageDataHandling<BodyTypeTuple>(), "Storage"));
+ //
+ //auto ccdID = forest->addBlockData(walberla::pe::ccd::createHashGridsDataHandling(globalBodyStorage, *storageId), "CCD");
+ //auto fcdID = forest->addBlockData(walberla::pe::fcd::createGenericFCDDataHandling<BodyTypeTuple, walberla::pe::fcd::AnalyticCollideFunctor>(), "FCD");
+
+ //cr = std::make_shared<walberla::pe::cr::HardContactSemiImplicitTimesteppingSolvers>(globalBodyStorage, forest, *storageId, ccdID, fcdID);
+ //cr->setMaxIterations(peParameter->maxPeIterations);
+ //cr->setRelaxationModel(walberla::pe::cr::HardContactSemiImplicitTimesteppingSolvers::ApproximateInelasticCoulombContactByDecoupling);
+ //cr->setRelaxationParameter(walberla::real_t(peParameter->relaxationParameter));
+ //cr->setGlobalLinearAcceleration(PeConverter::convert(peParameter->globalLinearAcceleration));
+
+ //for (auto blockIt = forest->begin(); blockIt != forest->end(); ++blockIt)
+ //{
+ // walberla::pe::ccd::ICCD* ccd = blockIt->getData< walberla::pe::ccd::ICCD >(ccdID);
+ // ccd->reloadBodies();
+ //}
+}
diff --git a/source/DemCoupling/physicsEngineAdapter/pe/PePhysicsEngineSolverAdapter.h b/source/DemCoupling/physicsEngineAdapter/pe/PePhysicsEngineSolverAdapter.h
index d4a5599fbce6f0ba1b72d01170f6e1f5def153e5..f79facea703a475b1d74ffacbc950f8aa9beb1c1 100644
--- a/source/DemCoupling/physicsEngineAdapter/pe/PePhysicsEngineSolverAdapter.h
+++ b/source/DemCoupling/physicsEngineAdapter/pe/PePhysicsEngineSolverAdapter.h
@@ -72,6 +72,8 @@ public:
walberla::pe::RigidBody* getPeGeoObject(walberla::id_t id);
void updateGeometry(std::shared_ptr<PhysicsEngineGeometryAdapter>) override;
std::shared_ptr< walberla::blockforest::BlockForest > getForest();
+ void saveToFile(const std::string& path);
+ void loadFromFile(const std::string& path);
private:
void initalizePeEnvironment();
diff --git a/source/VirtualFluidsBasic/basics/utilities/UbFileInputBinary.h b/source/VirtualFluidsBasic/basics/utilities/UbFileInputBinary.h
index 6b622c6a3a5de5730f307733175ea8b3c18aeb57..c410c822fbc97d5e7924efacf3987d8af657a0b6 100644
--- a/source/VirtualFluidsBasic/basics/utilities/UbFileInputBinary.h
+++ b/source/VirtualFluidsBasic/basics/utilities/UbFileInputBinary.h
@@ -64,6 +64,17 @@ public:
file.infile.read((char*)&data,sizeof(T));
return file;
}
+
+ template< typename T>
+ void readVector(std::vector<T>& v)
+ {
+ size_t size = v.size();
+ if (size > 0)
+ {
+ infile.read((char*)&v[0], sizeof(T)*size);
+ }
+ }
+
};
#endif
diff --git a/source/VirtualFluidsBasic/basics/utilities/UbFileOutputBinary.h b/source/VirtualFluidsBasic/basics/utilities/UbFileOutputBinary.h
index 13d94265c902abca8731382cfb08bd9e61d05e1c..76c50246793f61f294fe8305fb9fd46f0a807842 100644
--- a/source/VirtualFluidsBasic/basics/utilities/UbFileOutputBinary.h
+++ b/source/VirtualFluidsBasic/basics/utilities/UbFileOutputBinary.h
@@ -63,6 +63,16 @@ public:
file.outfile.write((char*)&data,sizeof(T));
return file;
}
+
+ template< typename T>
+ void writeVector(std::vector<T>& v)
+ {
+ size_t size = v.size();
+ if (size > 0)
+ {
+ outfile.write((char*)&v[0],sizeof(T)*size);
+ }
+ }
};
#endif
diff --git a/source/VirtualFluidsCore/BoundaryConditions/EqDensityBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/EqDensityBCAlgorithm.h
index d97af60907b7c97347a3d815ddf2cf0116e8a187..8c4f472bf1a8880bb4731006d9687dd13a2e3046 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/EqDensityBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/EqDensityBCAlgorithm.h
@@ -14,12 +14,5 @@ public:
SPtr<BCAlgorithm> clone();
void addDistributions(SPtr<DistributionArray3D> distributions);
void applyBC() override;
-private:
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // EqDensityBCAlgorithm_h__
diff --git a/source/VirtualFluidsCore/BoundaryConditions/HighViscosityNoSlipBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/HighViscosityNoSlipBCAlgorithm.h
index 264540fa2f95b0212dba4fe365aa990b0d3c13ba..aacc41a5a6aafb0bdd9312bc6ccf0e8fc618a874 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/HighViscosityNoSlipBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/HighViscosityNoSlipBCAlgorithm.h
@@ -13,15 +13,7 @@ public:
~HighViscosityNoSlipBCAlgorithm();
SPtr<BCAlgorithm> clone();
void addDistributions(SPtr<DistributionArray3D> distributions);
-
void applyBC() override;
-private:
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // HighViscosityNoSlipBCAlgorithm_h__
diff --git a/source/VirtualFluidsCore/BoundaryConditions/NoSlipBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/NoSlipBCAlgorithm.h
index 8f18a31f6454682ede88b125b35bf1f8e7c5ed8d..402e2e2bb70286533aafb7992a6186d28eb0e360 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/NoSlipBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/NoSlipBCAlgorithm.h
@@ -15,11 +15,5 @@ public:
void addDistributions(SPtr<DistributionArray3D> distributions);
void applyBC() override;
private:
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // NoSlipBCAlgorithm_h__
diff --git a/source/VirtualFluidsCore/BoundaryConditions/NonEqDensityBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/NonEqDensityBCAlgorithm.h
index 614c79bf34cb58f535d81af488205e762e6cc236..8cec91d6ec2e709dd23d0190c7422775b80319d0 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/NonEqDensityBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/NonEqDensityBCAlgorithm.h
@@ -14,12 +14,5 @@ public:
SPtr<BCAlgorithm> clone();
void addDistributions(SPtr<DistributionArray3D> distributions);
void applyBC() override;
-private:
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // NonEqDensityBCAlgorithm_h__
diff --git a/source/VirtualFluidsCore/BoundaryConditions/NonReflectingOutflowBCAlgorithm.cpp b/source/VirtualFluidsCore/BoundaryConditions/NonReflectingOutflowBCAlgorithm.cpp
index 31e80992a98ce0020d5019bed4d9db0362ca245f..149824f17c7d4472152571e617d9e766fe61a166 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/NonReflectingOutflowBCAlgorithm.cpp
+++ b/source/VirtualFluidsCore/BoundaryConditions/NonReflectingOutflowBCAlgorithm.cpp
@@ -9,7 +9,6 @@ NonReflectingOutflowBCAlgorithm::NonReflectingOutflowBCAlgorithm()
{
BCAlgorithm::type = BCAlgorithm::NonReflectingOutflowBCAlgorithm;
BCAlgorithm::preCollision = true;
- step=0;
}
//////////////////////////////////////////////////////////////////////////
NonReflectingOutflowBCAlgorithm::~NonReflectingOutflowBCAlgorithm()
diff --git a/source/VirtualFluidsCore/BoundaryConditions/NonReflectingOutflowBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/NonReflectingOutflowBCAlgorithm.h
index e2222f044c66c6ea1532fc4dbffb7f9aee3d1337..4effed280e073389eb1804983d6f7e10eb92aee0 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/NonReflectingOutflowBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/NonReflectingOutflowBCAlgorithm.h
@@ -14,13 +14,5 @@ public:
SPtr<BCAlgorithm> clone();
void addDistributions(SPtr<DistributionArray3D> distributions);
void applyBC() override;
-private:
- int step;
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // NonReflectingDensityBCAlgorithm_h__
diff --git a/source/VirtualFluidsCore/BoundaryConditions/SlipBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/SlipBCAlgorithm.h
index ee3f23a578ea50888662462f9eb83f1cf66d0da8..423cde915386649f63ac259bf76afc1b02381109 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/SlipBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/SlipBCAlgorithm.h
@@ -14,12 +14,5 @@ public:
SPtr<BCAlgorithm> clone();
void addDistributions(SPtr<DistributionArray3D> distributions);
void applyBC() override;
-private:
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // SlipBCAlgorithm_h__
diff --git a/source/VirtualFluidsCore/BoundaryConditions/ThinWallNoSlipBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/ThinWallNoSlipBCAlgorithm.h
index 68e0ff273fda0c0dbd99c7dfdbf8f90466e075b4..6f54bc7666ebd96bd9df5b3b07c698f29ff945fe 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/ThinWallNoSlipBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/ThinWallNoSlipBCAlgorithm.h
@@ -21,12 +21,5 @@ protected:
private:
int pass;
LBMReal fTemp[D3Q27System::ENDF + 1];
-
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // ThinWallNoSlipBCAlgorithm_h__
diff --git a/source/VirtualFluidsCore/BoundaryConditions/VelocityBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/VelocityBCAlgorithm.h
index dca549cd4e48467212c31e32d249009a61e353c7..c33ae68167a097eb1cdc298c7de3c7a074ff61ce 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/VelocityBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/VelocityBCAlgorithm.h
@@ -13,15 +13,7 @@ public:
~VelocityBCAlgorithm();
SPtr<BCAlgorithm> clone();
void addDistributions(SPtr<DistributionArray3D> distributions);
-
void applyBC() override;
-private:
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // VelocityBoundaryCondition_h__
diff --git a/source/VirtualFluidsCore/BoundaryConditions/VelocityWithDensityBCAlgorithm.h b/source/VirtualFluidsCore/BoundaryConditions/VelocityWithDensityBCAlgorithm.h
index 3de6b699071c6ca8014ff76f0de86def6283dc3d..25ec70608c1b8e54bc443207183b450e337b6884 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/VelocityWithDensityBCAlgorithm.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/VelocityWithDensityBCAlgorithm.h
@@ -20,12 +20,5 @@ public:
SPtr<BCAlgorithm> clone();
void addDistributions(SPtr<DistributionArray3D> distributions);
void applyBC();
-private:
- //friend class boost::serialization::access;
- //template<class Archive>
- //void serialize(Archive & ar, const unsigned int version)
- //{
- // ar & boost::serialization::base_object<BCAlgorithm>(*this);
- //}
};
#endif // NonReflectingVelocityBCAlgorithm_h__
diff --git a/source/VirtualFluidsCore/Grid/BasicCalculator.cpp b/source/VirtualFluidsCore/Grid/BasicCalculator.cpp
index b07fde93c2089420e102aabb1c679a60acce5ff2..8d7eafba0dca8905fe137398a604c09dc99e85bd 100644
--- a/source/VirtualFluidsCore/Grid/BasicCalculator.cpp
+++ b/source/VirtualFluidsCore/Grid/BasicCalculator.cpp
@@ -130,7 +130,8 @@ void BasicCalculator::calculate()
//now ghost nodes have actual values
}
UBLOG(logDEBUG1, "OMPCalculator::calculate() - stoped");
- }
catch (std::exception& e)
+ }
+ catch (std::exception& e)
{
UBLOG(logERROR, e.what());
UBLOG(logERROR, " step = "<<calcStep);
@@ -347,20 +348,20 @@ void BasicCalculator::applyPostCollisionBC(int startLevel, int maxInitLevel)
{
UBLOG(logERROR, e.what());
//UBLOG(logERROR, " step = "<<calcStep);
- throw;
- //exit(EXIT_FAILURE);
+ //throw;
+ exit(EXIT_FAILURE);
}
catch (std::string& s)
{
UBLOG(logERROR, s);
- throw;
- //exit(EXIT_FAILURE);
+ //throw;
+ exit(EXIT_FAILURE);
}
catch (...)
{
UBLOG(logERROR, "unknown exception");
- throw;
- //exit(EXIT_FAILURE);
+ //throw;
+ exit(EXIT_FAILURE);
}
}
diff --git a/source/VirtualFluidsCore/Parallel/Communicator.h b/source/VirtualFluidsCore/Parallel/Communicator.h
index fa7b98feb06b67e89b5469cb03af70ccd09df7b7..aa2e3d35d1371dbeb789e38d6a75848624b9df8e 100644
--- a/source/VirtualFluidsCore/Parallel/Communicator.h
+++ b/source/VirtualFluidsCore/Parallel/Communicator.h
@@ -41,9 +41,11 @@ public:
virtual void broadcast(int& value) = 0;
virtual void broadcast(float& value) = 0;
virtual void broadcast(double& value) = 0;
+ virtual void broadcast(long int& value) = 0;
virtual void broadcast(std::vector<int>& values) = 0;
virtual void broadcast(std::vector<float>& values) = 0;
virtual void broadcast(std::vector<double>& values) = 0;
+ virtual void broadcast(std::vector<long int>& values) = 0;
protected:
Communicator(){}
Communicator( const Communicator& ){}
diff --git a/source/VirtualFluidsCore/Parallel/MPICommunicator.cpp b/source/VirtualFluidsCore/Parallel/MPICommunicator.cpp
index f5c6090ba42427cfbbf25ee7381132b3cc518612..e2c32cdbcb8779150ea58e75eacd0bc710f5e925 100644
--- a/source/VirtualFluidsCore/Parallel/MPICommunicator.cpp
+++ b/source/VirtualFluidsCore/Parallel/MPICommunicator.cpp
@@ -208,6 +208,11 @@ void MPICommunicator::broadcast(std::vector<double>& values)
broadcast<double>(values);
}
//////////////////////////////////////////////////////////////////////////
+void MPICommunicator::broadcast(std::vector<long int>& values)
+{
+ broadcast<long int>(values);
+}
+//////////////////////////////////////////////////////////////////////////
void MPICommunicator::broadcast(int& value)
{
broadcast<int>(value);
@@ -223,5 +228,8 @@ void MPICommunicator::broadcast(double& value)
broadcast<double>(value);
}
//////////////////////////////////////////////////////////////////////////
-
+void MPICommunicator::broadcast(long int& value)
+{
+ broadcast<long int>(value);
+}
#endif
diff --git a/source/VirtualFluidsCore/Parallel/MPICommunicator.h b/source/VirtualFluidsCore/Parallel/MPICommunicator.h
index f21eda2db2f43763be35e627fb07eb28416b31e3..f57a7db4c9a1346239142517e54e28360d69cbf5 100644
--- a/source/VirtualFluidsCore/Parallel/MPICommunicator.h
+++ b/source/VirtualFluidsCore/Parallel/MPICommunicator.h
@@ -53,9 +53,11 @@ public:
void broadcast(int& value);
void broadcast(float& value);
void broadcast(double& value);
+ void broadcast(long int& value);
void broadcast(std::vector<int>& values);
void broadcast(std::vector<float>& values);
void broadcast(std::vector<double>& values);
+ void broadcast(std::vector<long int>& values);
template <class T>
std::vector<T> gather(std::vector<T>& values);
@@ -145,6 +147,7 @@ void MPICommunicator::broadcast(std::vector<T>& values)
if ((std::string)typeid(T).name()==(std::string)typeid(double).name()) mpiDataType = MPI_DOUBLE;
else if ((std::string)typeid(T).name()==(std::string)typeid(float).name()) mpiDataType = MPI_FLOAT;
else if ((std::string)typeid(T).name()==(std::string)typeid(int).name()) mpiDataType = MPI_INT;
+ else if ((std::string)typeid(T).name()==(std::string)typeid(long int).name()) mpiDataType = MPI_LONG_INT;
else throw UbException(UB_EXARGS, "no MpiDataType for T"+(std::string)typeid(T).name());
int rcount;
@@ -170,6 +173,7 @@ void MPICommunicator::broadcast(T& value)
if ((std::string)typeid(T).name() == (std::string)typeid(double).name()) mpiDataType = MPI_DOUBLE;
else if ((std::string)typeid(T).name() == (std::string)typeid(float).name()) mpiDataType = MPI_FLOAT;
else if ((std::string)typeid(T).name() == (std::string)typeid(int).name()) mpiDataType = MPI_INT;
+ else if ((std::string)typeid(T).name()==(std::string)typeid(long int).name()) mpiDataType = MPI_LONG_INT;
else throw UbException(UB_EXARGS, "no MpiDataType for T" + (std::string)typeid(T).name());
MPI_Bcast(&value, 1, mpiDataType, this->root, comm);
diff --git a/source/VirtualFluidsCore/Parallel/MetisPartitioner.cpp b/source/VirtualFluidsCore/Parallel/MetisPartitioner.cpp
index a234593cc7bb7a16e33cb331c58eecab5d910e42..31d23197e159933eb4dd6011853c1b45aaa4a38f 100644
--- a/source/VirtualFluidsCore/Parallel/MetisPartitioner.cpp
+++ b/source/VirtualFluidsCore/Parallel/MetisPartitioner.cpp
@@ -27,9 +27,9 @@ int MetisPartitioner::partition(int nofParts, MetisPartitioner::PartType ptype)
int rc;
idx_t nvtxs = (idx_t)xadj.size()-1; // number of nodes
idx_t ncon = (idx_t)vwgt.size()/nvtxs; // number Of node constraints;
-
part.resize(nvtxs);
- int edgecutCount = 0;
+ idx_t edgecutCount = 0;
+ idx_t nofPartsMetis = (idx_t)nofParts;
switch (ptype)
{
@@ -39,7 +39,7 @@ int MetisPartitioner::partition(int nofParts, MetisPartitioner::PartType ptype)
//else if( nofParts > 8 ) UBLOG(logWARNING, "MetisPartitioner::Recursive: !!!Warning!!! best for nofParts<=8 --> Kway is maybe a better option");
rc = METIS_PartGraphRecursive(&nvtxs, &ncon, &xadj[0], &adjncy[0],
- &vwgt[0], vsize, &adjwgt[0], &nofParts,
+ &vwgt[0], vsize, &adjwgt[0], &nofPartsMetis,
tpwgts, ubvec, options, &edgecutCount, &part[0]);
break;
case MetisPartitioner::KWAY:
@@ -48,7 +48,7 @@ int MetisPartitioner::partition(int nofParts, MetisPartitioner::PartType ptype)
//else if( nofParts < 9 ) UBLOG(logWARNING, "MetisPartitioner::Kway: !!!Warning!!! best for nofParts>8 --> Recursive is maybe a better option");
rc = METIS_PartGraphKway(&nvtxs, &ncon, &xadj[0], &adjncy[0],
- &vwgt[0], vsize, &adjwgt[0], &nofParts,
+ &vwgt[0], vsize, &adjwgt[0], &nofPartsMetis,
tpwgts, ubvec, options, &edgecutCount, &part[0]);
break;
}