From c1a8ef7c069a461f4830e4b25696d7570418545e Mon Sep 17 00:00:00 2001
From: Konstantin Kutscher <kutscher@irmb.tu-bs.de>
Date: Thu, 5 Jul 2018 18:36:36 +0200
Subject: [PATCH] -fixes bad allocation issue in CreateDemObjectsCoProcessor
-fixes destructor crash in WriteDemObjectsCoProcessor
---
source/Applications/Thermoplast/config.txt | 25 +-
.../Applications/Thermoplast/thermoplast.cpp | 539 ++++++++----------
.../CreateDemObjectsCoProcessor.cpp | 31 +-
source/DemCoupling/DemCoProcessor.cpp | 2 +-
.../RestartDemObjectsCoProcessor.cpp | 23 +-
.../WriteDemObjectsCoProcessor.cpp | 2 +-
.../DemCoupling/WriteDemObjectsCoProcessor.h | 4 +-
7 files changed, 278 insertions(+), 348 deletions(-)
diff --git a/source/Applications/Thermoplast/config.txt b/source/Applications/Thermoplast/config.txt
index d8e5a35d2..be02fd59a 100644
--- a/source/Applications/Thermoplast/config.txt
+++ b/source/Applications/Thermoplast/config.txt
@@ -1,24 +1,31 @@
#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
+boundingBox = 60 1370 10 190 1530 320 #test bb
+
+#boundingBox = 60 0 10 190 1530 2320 #production bb
+
blocknx = 10 10 10
#blocknx = 300 420 320
-endTime = 100
-outTime = 100
+endTime = 10
+outTime = 10
availMem = 25e9
uLB = 0.1
#PE parameters
-peMinOffset = 0 0 0 #46 2 2
-peMaxOffset = 0 0 0 #-8 -25 -2
-sphereTime = 100
+peMinOffset = 46 2 2
+peMaxOffset = -8 -25 -2
+sphereTime = 10
#geometry files
pathGeo = d:/Projects/ThermoPlast/SimGeo
michel = /michel.stl
plexiglas = /plexiglas.stl
-pathOut = g:/temp/thermoplastCluster3
\ No newline at end of file
+pathOut = g:/temp/thermoplastCluster3
+
+
+#restart
+cpStart = 200
+cpStep = 200
+restart = false
\ No newline at end of file
diff --git a/source/Applications/Thermoplast/thermoplast.cpp b/source/Applications/Thermoplast/thermoplast.cpp
index 275e5b877..c10d50ef3 100644
--- a/source/Applications/Thermoplast/thermoplast.cpp
+++ b/source/Applications/Thermoplast/thermoplast.cpp
@@ -31,7 +31,7 @@
#include <WriteDemObjectsCoProcessor.h>
#include "CreateDemObjectsCoProcessor.h"
-#include "RestartDemCouplingCoProcessor.h"
+#include "RestartDemObjectsCoProcessor.h"
using namespace std;
@@ -85,69 +85,149 @@ std::shared_ptr<DemCoProcessor> makePeCoProcessor(SPtr<Grid3D> grid, SPtr<Commun
return std::make_shared<DemCoProcessor>(grid, peScheduler, comm, forceCalculator, peSolver);
}
-void createGridObjects()
-{
-
-}
-
//pipe flow with forcing
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");
+
+ 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 michel = config.getValue<string>("michel");
+ string plexiglas = config.getValue<string>("plexiglas");
+ double sphereTime = config.getValue<double>("sphereTime");
+
+ double cpStart = config.getValue<double>("cpStart");
+ double cpStep = config.getValue<double>("cpStep");
+ bool restart = config.getValue<bool>("restart");
+
+ peMinOffset = config.getVector<double>("peMinOffset");
+ peMaxOffset = config.getVector<double>("peMaxOffset");
+
+ pathOut = config.getValue<string>("pathOut");
+ pathGeo = config.getValue<string>("pathGeo");
+
+ //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;
+
+ //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);
+
+ //boundary conditions definition
+ //////////////////////////////////////////////////////////////////////////////
+ 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 BC
+ 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()));
+
+ //boundary conditions visitor
+ SPtr<BoundaryConditionsBlockVisitor> bcVisitor(new BoundaryConditionsBlockVisitor());
+ bcVisitor->addBC(noSlipBCAdapter);
+ bcVisitor->addBC(inflowAdapter);
+ bcVisitor->addBC(outflowAdapter);
+ bcVisitor->addBC(velocityBcParticleAdapter);
+ //////////////////////////////////////////////////////////////////////////////////
+
+ //LBM kernel definition
+ SPtr<LBMKernel> kernel;
+ kernel = SPtr<LBMKernel>(new IncompressibleCumulantLBMKernel());
+ SPtr<BCProcessor> bcProc(new BCProcessor());
+ kernel->setBCProcessor(bcProc);
+
+ //blocks generating
+ 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());
+ GenBlocksGridVisitor genBlocks(gridCube);
+ grid->accept(genBlocks);
+
+ //PE initialization
+ 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->setBlockVisitor(bcVisitor);
+
+ if (myid == 0)
{
- 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());
+ 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");
+ }
+ 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());
+
+ if (!restart)
+ {
//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());
@@ -164,60 +244,7 @@ void thermoplast(string configname)
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());
@@ -231,23 +258,9 @@ void thermoplast(string configname)
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));
@@ -264,51 +277,15 @@ void thermoplast(string configname)
//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(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
@@ -341,40 +318,18 @@ void thermoplast(string configname)
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));
@@ -386,83 +341,113 @@ void thermoplast(string configname)
}
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));
+ ////////////////////////////////////////////////////////////////////////////
+ ////generating spheres
+ //UBLOG(logINFO, "generating spheres - start, rank="<<myid);
+ SPtr<UbScheduler> sphereScheduler(new UbScheduler(sphereTime));
+ SPtr<CreateDemObjectsCoProcessor> createSphereCoProcessor(new CreateDemObjectsCoProcessor(grid, sphereScheduler, demCoProcessor, sphereMaterial));
+ //UBLOG(logINFO, "generating spheres - stop, rank="<<myid);
+
+ //restart
+ //UBLOG(logINFO, "restart definition - start, rank="<<myid);
+ SPtr<UbScheduler> restartSch(new UbScheduler(cpStep, cpStart));
+ SPtr<MPIIORestartCoProcessor> restartCoProcessor(new MPIIORestartCoProcessor(grid, restartSch, pathOut, comm));
+ restartCoProcessor->setLBMKernel(kernel);
+ restartCoProcessor->setBCProcessor(bcProc);
+ SPtr<RestartDemObjectsCoProcessor> restartDemObjectsCoProcessor(new RestartDemObjectsCoProcessor(grid, restartSch, pathOut, demCoProcessor, createSphereCoProcessor, radius, comm));
+ //UBLOG(logINFO, "restart definition - stop, rank="<<myid);
+
+ if (restart)
+ {
+ int startStep = restartCoProcessor->readCpTimeStep();
+ restartCoProcessor->restart(startStep);
+ restartDemObjectsCoProcessor->restart(startStep);
+ }
- SPtr<WriteDemObjectsCoProcessor> writeDemObjectsCoProcessor(new WriteDemObjectsCoProcessor(grid, visSch, pathOut, SPtr<WbWriter>(WbWriterVtkXmlBinary::getInstance()), demCoProcessor, comm));
- writeDemObjectsCoProcessor->process(0);
+ //set connectors
+ //UBLOG(logINFO, "set connectors - start, rank="<<myid);
+ InterpolationProcessorPtr iProcessor(new IncompressibleOffsetInterpolationProcessor());
+ SetConnectorsBlockVisitor setConnsVisitor(comm, true, D3Q27System::ENDDIR, nuLB, iProcessor);
+ grid->accept(setConnsVisitor);
+ //UBLOG(logINFO, "set connectors - stop, rank="<<myid);
+
+ //BC visitor
+ //UBLOG(logINFO, "BC visitor - start, rank="<<myid);
+ grid->accept(*bcVisitor.get());
+ //UBLOG(logINFO, "BC visitor - stop, rank="<<myid);
+
+ //sphere prototypes
+ //UBLOG(logINFO, "sphere prototypes - start, rank="<<myid);
+ double d = 2.0*radius;
+ Vector3D origin(g_minX1+peMinOffset[0]+radius, geoInflow1->getX2Minimum()+4.0*d, geoInflow1->getX3Minimum()+2.0*d);
+ 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));
+ }
+ //UBLOG(logINFO, "sphere prototypes - stop, rank="<<myid);
+
+ //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, 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));
+
+ //start simulation
+ omp_set_num_threads(numOfThreads);
+ SPtr<UbScheduler> stepGhostLayer(peScheduler);
+ SPtr<Calculator> calculator(new BasicCalculator(grid, stepGhostLayer, 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");
+}
- //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");
+//////////////////////////////////////////////////////////////////////////
+int main(int argc, char* argv[])
+{
+ try
+ {
+ //Sleep(30000);
+ walberla::Environment env(argc, argv);
+ if (argv!=NULL)
+ {
+ if (argv[1]!=NULL)
+ {
+ thermoplast(string(argv[1]));
+ }
+ else
+ {
+ cout<<"Configuration file must be set!: "<<argv[0]<<" <config file>"<<endl<<std::flush;
+ }
+ }
+ return 0;
}
catch (std::exception& e)
{
@@ -477,39 +462,3 @@ void thermoplast(string configname)
UBLOG(logERROR, "unknown exception");
}
}
-
-
-//////////////////////////////////////////////////////////////////////////
-int main(int argc, char* argv[])
-{
- //try
- //{
- //Sleep(30000);
- walberla::Environment env(argc, argv);
-
- if (argv!=NULL)
- {
- if (argv[1]!=NULL)
- {
- thermoplast(string(argv[1]));
- }
- 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");
- //}
-}
diff --git a/source/DemCoupling/CreateDemObjectsCoProcessor.cpp b/source/DemCoupling/CreateDemObjectsCoProcessor.cpp
index cab9e39d7..f6b2c13ca 100644
--- a/source/DemCoupling/CreateDemObjectsCoProcessor.cpp
+++ b/source/DemCoupling/CreateDemObjectsCoProcessor.cpp
@@ -20,8 +20,7 @@
CreateDemObjectsCoProcessor::CreateDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, SPtr<DemCoProcessor> demCoProcessor, SPtr<PhysicsEngineMaterialAdapter> demObjectMaterial) :
CoProcessor(grid, s),
demCoProcessor(demCoProcessor),
- demObjectMaterial(demObjectMaterial),
- initalVelocity(initalVelocity)
+ demObjectMaterial(demObjectMaterial)
{
mu::Parser fct;
fct.SetExpr("U");
@@ -40,34 +39,6 @@ void CreateDemObjectsCoProcessor::process(double step)
if (scheduler->isDue(step))
{
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));
-
- //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);
- //}
}
}
//////////////////////////////////////////////////////////////////////////
diff --git a/source/DemCoupling/DemCoProcessor.cpp b/source/DemCoupling/DemCoProcessor.cpp
index ce4de51cd..34686abff 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();
}
diff --git a/source/DemCoupling/RestartDemObjectsCoProcessor.cpp b/source/DemCoupling/RestartDemObjectsCoProcessor.cpp
index 2e04192a4..da7e0c24c 100644
--- a/source/DemCoupling/RestartDemObjectsCoProcessor.cpp
+++ b/source/DemCoupling/RestartDemObjectsCoProcessor.cpp
@@ -25,27 +25,29 @@ void RestartDemObjectsCoProcessor::process(double step)
{
int istep = static_cast<int>(step);
- write(istep);
+ if (comm->isRoot())
+ UBLOG(logINFO, "RestartDemObjectsCoProcessor::write step: " << istep);
- if (comm->isRoot())
- UBLOG(logINFO, "RestartDemObjectsCoProcessor write step: " << istep);
+ write(istep);
}
}
void RestartDemObjectsCoProcessor::restart(double step)
{
+ if (comm->isRoot())
+ UBLOG(logINFO, "RestartDemObjectsCoProcessor::read step: " << (int)step);
+
read((int)step);
}
void RestartDemObjectsCoProcessor::write(int step)
{
- if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor write:start ");
+ 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());
-
+ //TODO implement getherv
std::vector<double> rvalues;
comm->allGather(p, rvalues);
@@ -56,13 +58,14 @@ void RestartDemObjectsCoProcessor::write(int step)
UbFileOutputBinary fo(filePath);
fo.writeInteger((int)rvalues.size());
fo.writeVector<double>(rvalues);
- UBLOG(logINFO, "RestartDemObjectsCoProcessor size: " << rvalues.size());
+ UBLOG(logINFO, "RestartDemObjectsCoProcessor: number of objects = " << rvalues.size()/6);
}
- if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor write:stop ");
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor::write stop ");
}
void RestartDemObjectsCoProcessor::read(int step)
{
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor::read start ");
std::vector<double> p;
if (comm->isRoot())
@@ -73,11 +76,10 @@ void RestartDemObjectsCoProcessor::read(int step)
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());
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor::read number of objects = " << p.size()/6);
createDemObjectsCoProcessor->clearGeoObjects();
@@ -93,4 +95,5 @@ void RestartDemObjectsCoProcessor::read(int step)
createDemObjectsCoProcessor->clearGeoObjects();
+ if (comm->isRoot()) UBLOG(logINFO, "RestartDemObjectsCoProcessor::read stop ");
}
diff --git a/source/DemCoupling/WriteDemObjectsCoProcessor.cpp b/source/DemCoupling/WriteDemObjectsCoProcessor.cpp
index 682433a2e..130cd3027 100644
--- a/source/DemCoupling/WriteDemObjectsCoProcessor.cpp
+++ b/source/DemCoupling/WriteDemObjectsCoProcessor.cpp
@@ -14,7 +14,7 @@ WriteDemObjectsCoProcessor::WriteDemObjectsCoProcessor()
}
//////////////////////////////////////////////////////////////////////////
-WriteDemObjectsCoProcessor::WriteDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, const std::string& path, SPtr<WbWriter> writer, SPtr<DemCoProcessor> demCoProcessor, SPtr<Communicator> comm)
+WriteDemObjectsCoProcessor::WriteDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, const std::string& path, WbWriter* const writer, SPtr<DemCoProcessor> demCoProcessor, SPtr<Communicator> comm)
: CoProcessor(grid, s),
path(path),
writer(writer),
diff --git a/source/DemCoupling/WriteDemObjectsCoProcessor.h b/source/DemCoupling/WriteDemObjectsCoProcessor.h
index 81b00b745..93d87bcf8 100644
--- a/source/DemCoupling/WriteDemObjectsCoProcessor.h
+++ b/source/DemCoupling/WriteDemObjectsCoProcessor.h
@@ -21,13 +21,13 @@ class WriteDemObjectsCoProcessor : public CoProcessor
{
public:
WriteDemObjectsCoProcessor();
- WriteDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, const std::string& path, SPtr<WbWriter> const writer, SPtr<DemCoProcessor> demCoProcessor, SPtr<Communicator> comm);
+ WriteDemObjectsCoProcessor(SPtr<Grid3D> grid, SPtr<UbScheduler> s, const std::string& path, WbWriter* const writer, SPtr<DemCoProcessor> demCoProcessor, SPtr<Communicator> comm);
~WriteDemObjectsCoProcessor() {}
void process(double step) override;
private:
std::string path;
- SPtr<WbWriter> writer;
+ WbWriter* writer;
SPtr<Communicator> comm;
SPtr<DemCoProcessor> demCoProcessor;
};
--
GitLab