diff --git a/apps/cpu/GyroidsRow/GyroidsRow.cpp b/apps/cpu/GyroidsRow/GyroidsRow.cpp index c6a622d9761ca6716fc9d0b50f8d5c247f1a9897..d9d16743cc15598134b4296490c26e27b7460aa6 100644 --- a/apps/cpu/GyroidsRow/GyroidsRow.cpp +++ b/apps/cpu/GyroidsRow/GyroidsRow.cpp @@ -42,318 +42,311 @@ using namespace std; using namespace vf::lbm::dir; using namespace vf::basics::constant; -void run(string configname) +void run(const vf::basics::ConfigurationFile& config) { - try { - vf::basics::ConfigurationFile config; - config.load(configname); - - string pathname = config.getValue<string>("pathname"); - int numOfThreads = config.getValue<int>("numOfThreads"); - vector<int> blocknx = config.getVector<int>("blocknx"); - double endTime = config.getValue<double>("endTime"); - double outTime = config.getValue<double>("outTime"); - double availMem = config.getValue<double>("availMem"); - double nu = config.getValue<double>("nu"); - double dx = config.getValue<double>("dx"); - double UnitEdgeLength = config.getValue<double>("UnitEdgeLength"); - double Re = config.getValue<double>("Re"); - double vx = config.getValue<double>("vx"); - vector<double> length = config.getVector<double>("length"); - //double timeAvStart = config.getValue<double>("timeAvStart"); - //double timeAvStop = config.getValue<double>("timeAvStop"); - vector<double> TPMSL = config.getVector<double>("TPMSL"); - vector<double> TPMSOrigin = config.getVector<double>("TPMSOrigin"); - vector<double> gridCubeOrigin = config.getVector<double>("gridCubeOrigin"); - int refineLevel = config.getValue<int>("refineLevel"); - bool logToFile = config.getValue<bool>("logToFile"); - double restartStep = config.getValue<double>("restartStep"); - double cpStart = config.getValue<double>("cpStart"); - double cpStep = config.getValue<double>("cpStep"); - bool newStart = config.getValue<bool>("newStart"); - - SPtr<vf::parallel::Communicator> comm = vf::parallel::MPICommunicator::getInstance(); - int myid = comm->getProcessID(); - - if (logToFile) { + string pathname = config.getValue<string>("pathname"); + int numOfThreads = config.getValue<int>("numOfThreads"); + vector<int> blocknx = config.getVector<int>("blocknx"); + double endTime = config.getValue<double>("endTime"); + double outTime = config.getValue<double>("outTime"); + double availMem = config.getValue<double>("availMem"); + double nu = config.getValue<double>("nu"); + double dx = config.getValue<double>("dx"); + double UnitEdgeLength = config.getValue<double>("UnitEdgeLength"); + double Re = config.getValue<double>("Re"); + double vx = config.getValue<double>("vx"); + vector<double> length = config.getVector<double>("length"); + //double timeAvStart = config.getValue<double>("timeAvStart"); + //double timeAvStop = config.getValue<double>("timeAvStop"); + vector<double> TPMSL = config.getVector<double>("TPMSL"); + vector<double> TPMSOrigin = config.getVector<double>("TPMSOrigin"); + vector<double> gridCubeOrigin = config.getVector<double>("gridCubeOrigin"); + int refineLevel = config.getValue<int>("refineLevel"); + bool logToFile = config.getValue<bool>("logToFile"); + double restartStep = config.getValue<double>("restartStep"); + double cpStart = config.getValue<double>("cpStart"); + double cpStep = config.getValue<double>("cpStep"); + bool newStart = config.getValue<bool>("newStart"); + + SPtr<vf::parallel::Communicator> comm = vf::parallel::MPICommunicator::getInstance(); + int myid = comm->getProcessID(); + + if (logToFile) { #if defined(__unix__) - if (myid == 0) { - const char *str = pathname.c_str(); - mkdir(str, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH); - } + if (myid == 0) { + const char *str = pathname.c_str(); + mkdir(str, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH); + } #endif - if (myid == 0) { - stringstream logFilename; - logFilename << pathname + "/logfile" + UbSystem::toString(UbSystem::getTimeStamp()) + ".txt"; - UbLog::output_policy::setStream(logFilename.str()); - } + if (myid == 0) { + stringstream logFilename; + logFilename << pathname + "/logfile" + UbSystem::toString(UbSystem::getTimeStamp()) + ".txt"; + UbLog::output_policy::setStream(logFilename.str()); } + } - SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter()); - - - - //////////////////////////////////////////////////////////////////////// - // BC Adapter - SPtr<BC> tpmsNoslipAdapter(new NoSlipBC()); - SPtr<BC> xMinApr(new VelocityBC(vx, 0., BCFunction::INFCONST, 0., 0., BCFunction::INFCONST, 0.,0., BCFunction::INFCONST)); - SPtr<BC> xMaxApr(new PressureBC(0.)); - SPtr<BC> zMinApr(new NoSlipBC()); - SPtr<BC> zMaxApr(new NoSlipBC()); - - tpmsNoslipAdapter->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); - xMinApr->setBCStrategy(SPtr<BCStrategy>(new VelocityNonReflecting(c1o2))); - xMaxApr->setBCStrategy(SPtr<BCStrategy>(new OutflowNonReflectingWithPressure(c1o100))); - zMinApr->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); - zMaxApr->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); + SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter()); + + + + //////////////////////////////////////////////////////////////////////// + // BC Adapter + SPtr<BC> tpmsNoslipAdapter(new NoSlipBC()); + SPtr<BC> xMinApr(new VelocityBC(vx, 0., BCFunction::INFCONST, 0., 0., BCFunction::INFCONST, 0.,0., BCFunction::INFCONST)); + SPtr<BC> xMaxApr(new PressureBC(0.)); + SPtr<BC> zMinApr(new NoSlipBC()); + SPtr<BC> zMaxApr(new NoSlipBC()); + + tpmsNoslipAdapter->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); + xMinApr->setBCStrategy(SPtr<BCStrategy>(new VelocityNonReflecting(c1o2))); + xMaxApr->setBCStrategy(SPtr<BCStrategy>(new OutflowNonReflectingWithPressure(c1o100))); + zMinApr->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); + zMaxApr->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); + + ////////////////////////////////////////////////////////////////////////234 + // BC visitor + BoundaryConditionsBlockVisitor bcVisitor; + //////////////////////////////////////////////////////////////////////// + // grid, kernel and BCProcessor + SPtr<Grid3D> grid(new Grid3D(comm)); + SPtr<LBMKernel> kernel; + kernel = SPtr<LBMKernel>(new K15CompressibleNavierStokes());; + SPtr<BCSet> bcProc(new BCSet()); + kernel->setBCSet(bcProc); + SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelIntersected, d00M, MetisPartitioner::RECURSIVE)); + + ////////////////////////////////////////////////////////////////////////// + // restart + SPtr<UbScheduler> mSch(new UbScheduler(cpStep, cpStart)); + SPtr<MPIIOMigrationSimulationObserver> migSimulationObserver( + new MPIIOMigrationSimulationObserver(grid, mSch,metisVisitor, pathname + "/mig", comm)); + migSimulationObserver->setLBMKernel(kernel); + migSimulationObserver->setBCSet(bcProc); + ////////////////////////////////////////////////////////////////////////// + + if (newStart) { + //GbGyroidThirdOrderPtr tpms; + // tpms = GbGyroidThirdOrderPtr(new GbGyroidThirdOrder(TPMSOrigin[0], TPMSOrigin[1], TPMSOrigin[2], + // TPMSOrigin[0] + TPMSL[0], + // TPMSOrigin[1] + TPMSL[1], + // TPMSOrigin[2] + TPMSL[2], + // UnitEdgeLength, dx, 2.5e-4)); + GbGyroidThirdOrderLongPtr tpms; + tpms = GbGyroidThirdOrderLongPtr(new GbGyroidThirdOrderLong(TPMSOrigin[0], TPMSOrigin[1], TPMSOrigin[2], + TPMSOrigin[0] + TPMSL[0], + TPMSOrigin[1] + TPMSL[1], + TPMSOrigin[2] + TPMSL[2], + UnitEdgeLength, dx, 2.5e-4)); + + double g_minX1 = gridCubeOrigin[0]; + double g_minX2 = gridCubeOrigin[1]; + double g_minX3 = gridCubeOrigin[2]; + + double g_maxX1 = gridCubeOrigin[0] + length[0]; + double g_maxX2 = gridCubeOrigin[1] + length[1]; + double g_maxX3 = gridCubeOrigin[2] + length[2]; + + SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3)); + if (myid == 0) + GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", + WbWriterVtkXmlBinary::getInstance()); - ////////////////////////////////////////////////////////////////////////234 - // BC visitor - BoundaryConditionsBlockVisitor bcVisitor; - //////////////////////////////////////////////////////////////////////// - // grid, kernel and BCProcessor - SPtr<Grid3D> grid(new Grid3D(comm)); - SPtr<LBMKernel> kernel; - kernel = SPtr<LBMKernel>(new K15CompressibleNavierStokes());; - SPtr<BCSet> bcProc(new BCSet()); - kernel->setBCSet(bcProc); - SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelIntersected, d00M, MetisPartitioner::RECURSIVE)); + + SPtr<GbCuboid3D> spongecube(new GbCuboid3D(TPMSOrigin[0] + TPMSL[0], g_minX2 - dx, g_minX3 - dx, + g_maxX1 + dx, g_maxX2 + dx, g_maxX3 + dx)); + if (myid == 0) + GbSystem3D::writeGeoObject(spongecube.get(), pathname + "/geo/spongecube", + WbWriterVtkXmlBinary::getInstance()); + if (myid == 0) { + // UBLOG(logINFO,"rho = " << rhoLB ); + UBLOG(logINFO, "nu = " << nu); + UBLOG(logINFO, "Re = " << Re); + UBLOG(logINFO, "vx = " << vx); + UBLOG(logINFO, "dx = " << dx); + UBLOG(logINFO, "Preprocess - start"); + } - ////////////////////////////////////////////////////////////////////////// - // restart - SPtr<UbScheduler> mSch(new UbScheduler(cpStep, cpStart)); - SPtr<MPIIOMigrationSimulationObserver> migSimulationObserver( - new MPIIOMigrationSimulationObserver(grid, mSch,metisVisitor, pathname + "/mig", comm)); - migSimulationObserver->setLBMKernel(kernel); - migSimulationObserver->setBCSet(bcProc); - ////////////////////////////////////////////////////////////////////////// + grid->setDeltaX(dx); + grid->setBlockNX(blocknx[0], blocknx[1], blocknx[2]); + grid->setPeriodicX1(false); + grid->setPeriodicX2(true); + grid->setPeriodicX3(false); - if (newStart) { - //GbGyroidThirdOrderPtr tpms; - // tpms = GbGyroidThirdOrderPtr(new GbGyroidThirdOrder(TPMSOrigin[0], TPMSOrigin[1], TPMSOrigin[2], - // TPMSOrigin[0] + TPMSL[0], - // TPMSOrigin[1] + TPMSL[1], - // TPMSOrigin[2] + TPMSL[2], - // UnitEdgeLength, dx, 2.5e-4)); - GbGyroidThirdOrderLongPtr tpms; - tpms = GbGyroidThirdOrderLongPtr(new GbGyroidThirdOrderLong(TPMSOrigin[0], TPMSOrigin[1], TPMSOrigin[2], - TPMSOrigin[0] + TPMSL[0], - TPMSOrigin[1] + TPMSL[1], - TPMSOrigin[2] + TPMSL[2], - UnitEdgeLength, dx, 2.5e-4)); - - double g_minX1 = gridCubeOrigin[0]; - double g_minX2 = gridCubeOrigin[1]; - double g_minX3 = gridCubeOrigin[2]; - - double g_maxX1 = gridCubeOrigin[0] + length[0]; - double g_maxX2 = gridCubeOrigin[1] + length[1]; - double g_maxX3 = gridCubeOrigin[2] + length[2]; - - SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3)); - if (myid == 0) - GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", - WbWriterVtkXmlBinary::getInstance()); - - - SPtr<GbCuboid3D> spongecube(new GbCuboid3D(TPMSOrigin[0] + TPMSL[0], g_minX2 - dx, g_minX3 - dx, - g_maxX1 + dx, g_maxX2 + dx, g_maxX3 + dx)); - if (myid == 0) - GbSystem3D::writeGeoObject(spongecube.get(), pathname + "/geo/spongecube", - WbWriterVtkXmlBinary::getInstance()); - if (myid == 0) { - // UBLOG(logINFO,"rho = " << rhoLB ); - UBLOG(logINFO, "nu = " << nu); - UBLOG(logINFO, "Re = " << Re); - UBLOG(logINFO, "vx = " << vx); - UBLOG(logINFO, "dx = " << dx); - UBLOG(logINFO, "Preprocess - start"); - } + GenBlocksGridVisitor genBlocks(gridCube); + grid->accept(genBlocks); - grid->setDeltaX(dx); - grid->setBlockNX(blocknx[0], blocknx[1], blocknx[2]); - grid->setPeriodicX1(false); - grid->setPeriodicX2(true); - grid->setPeriodicX3(false); + SPtr<SimulationObserver> ppblocks(new WriteBlocksSimulationObserver(grid, SPtr<UbScheduler>(new UbScheduler(1)), pathname, + WbWriterVtkXmlBinary::getInstance(), comm)); - GenBlocksGridVisitor genBlocks(gridCube); - grid->accept(genBlocks); + ppblocks->update(0); - SPtr<SimulationObserver> ppblocks(new WriteBlocksSimulationObserver(grid, SPtr<UbScheduler>(new UbScheduler(1)), pathname, - WbWriterVtkXmlBinary::getInstance(), comm)); - ppblocks->update(0); + GbCuboid3DPtr xMin( new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_minX1, g_maxX2 + dx, g_maxX3 + dx)); + GbCuboid3DPtr xMax(new GbCuboid3D(g_maxX1 , g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_maxX3 + dx)); + GbCuboid3DPtr zMin(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_minX3)); + GbCuboid3DPtr zMax(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_maxX3, g_maxX1 + dx, g_maxX2 + dx, g_maxX3 + dx)); + if (myid == 0) + GbSystem3D::writeGeoObject(xMin.get(), pathname + "/geo/xMin", WbWriterVtkXmlBinary::getInstance()); + if (myid == 0) + GbSystem3D::writeGeoObject(xMax.get(), pathname + "/geo/xMax", WbWriterVtkXmlBinary::getInstance()); + if (myid == 0) + GbSystem3D::writeGeoObject(zMin.get(), pathname + "/geo/zMin", WbWriterVtkXmlBinary::getInstance()); + if (myid == 0) + GbSystem3D::writeGeoObject(zMax.get(), pathname + "/geo/zMax", WbWriterVtkXmlBinary::getInstance()); - GbCuboid3DPtr xMin( new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_minX1, g_maxX2 + dx, g_maxX3 + dx)); - GbCuboid3DPtr xMax(new GbCuboid3D(g_maxX1 , g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_maxX3 + dx)); - GbCuboid3DPtr zMin(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_minX3)); - GbCuboid3DPtr zMax(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_maxX3, g_maxX1 + dx, g_maxX2 + dx, g_maxX3 + dx)); + SPtr<D3Q27Interactor> tpmsInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(tpms, grid, tpmsNoslipAdapter, Interactor3D::SOLID, Interactor3D::POINTS)); + SPtr<D3Q27Interactor> xMinInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(xMin, grid, xMinApr, Interactor3D::SOLID, Interactor3D::POINTS)); + SPtr<D3Q27Interactor> xMaxInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(xMax, grid, xMaxApr, Interactor3D::SOLID, Interactor3D::POINTS)); + SPtr<D3Q27Interactor> zMinInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(zMin, grid, zMinApr, Interactor3D::SOLID, Interactor3D::POINTS)); + SPtr<D3Q27Interactor> zMaxInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(zMax, grid, zMaxApr, Interactor3D::SOLID, Interactor3D::POINTS)); - if (myid == 0) - GbSystem3D::writeGeoObject(xMin.get(), pathname + "/geo/xMin", WbWriterVtkXmlBinary::getInstance()); - if (myid == 0) - GbSystem3D::writeGeoObject(xMax.get(), pathname + "/geo/xMax", WbWriterVtkXmlBinary::getInstance()); - if (myid == 0) - GbSystem3D::writeGeoObject(zMin.get(), pathname + "/geo/zMin", WbWriterVtkXmlBinary::getInstance()); - if (myid == 0) - GbSystem3D::writeGeoObject(zMax.get(), pathname + "/geo/zMax", WbWriterVtkXmlBinary::getInstance()); + InteractorsHelper intHelper(grid, metisVisitor,false); - SPtr<D3Q27Interactor> tpmsInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(tpms, grid, tpmsNoslipAdapter, Interactor3D::SOLID, Interactor3D::POINTS)); - SPtr<D3Q27Interactor> xMinInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(xMin, grid, xMinApr, Interactor3D::SOLID, Interactor3D::POINTS)); - SPtr<D3Q27Interactor> xMaxInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(xMax, grid, xMaxApr, Interactor3D::SOLID, Interactor3D::POINTS)); - SPtr<D3Q27Interactor> zMinInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(zMin, grid, zMinApr, Interactor3D::SOLID, Interactor3D::POINTS)); - SPtr<D3Q27Interactor> zMaxInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(zMax, grid, zMaxApr, Interactor3D::SOLID, Interactor3D::POINTS)); + intHelper.addInteractor(tpmsInt); + intHelper.addInteractor(zMinInt); + intHelper.addInteractor(zMaxInt); - InteractorsHelper intHelper(grid, metisVisitor,false); + intHelper.addInteractor(xMinInt); + intHelper.addInteractor(xMaxInt); - intHelper.addInteractor(tpmsInt); - intHelper.addInteractor(zMinInt); - intHelper.addInteractor(zMaxInt); - intHelper.addInteractor(xMinInt); - intHelper.addInteractor(xMaxInt); + intHelper.selectBlocks(); + + // domain decomposition for threads + //PQueuePartitioningGridVisitor pqPartVisitor(numOfThreads); + //grid->accept(pqPartVisitor); + ppblocks->update(0); + ppblocks.reset(); - intHelper.selectBlocks(); - - // domain decomposition for threads - //PQueuePartitioningGridVisitor pqPartVisitor(numOfThreads); - //grid->accept(pqPartVisitor); - - ppblocks->update(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++) { - 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"); + ////////////////////////////////////////////////////////////////////////// + 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++) { + 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"); + } + ////////////////////////////////////////////////////////////////////////// - SetKernelBlockVisitor kernelVisitor(kernel, nu); - grid->accept(kernelVisitor); + SetKernelBlockVisitor kernelVisitor(kernel, nu); + grid->accept(kernelVisitor); - intHelper.setBC(); + intHelper.setBC(); - SpongeLayerBlockVisitor spongeLayerVisitor(spongecube, kernel, nu, dP00); - grid->accept(spongeLayerVisitor); + SpongeLayerBlockVisitor spongeLayerVisitor(spongecube, kernel, nu, dP00); + grid->accept(spongeLayerVisitor); - grid->accept(bcVisitor); + grid->accept(bcVisitor); - // initialization of distributions - InitDistributionsBlockVisitor initVisitor; - grid->accept(initVisitor); + // initialization of distributions + InitDistributionsBlockVisitor initVisitor; + grid->accept(initVisitor); - // boundary conditions grid - { - SPtr<UbScheduler> geoSch(new UbScheduler(1)); - SPtr<SimulationObserver> ppgeo(new WriteBoundaryConditionsSimulationObserver(grid, geoSch, pathname, WbWriterVtkXmlBinary::getInstance(), comm)); - ppgeo->update(0); - ppgeo.reset(); - } - if (myid == 0) - UBLOG(logINFO, "Preprocess - end"); - } - else + // boundary conditions grid { - if (myid == 0) { - UBLOG(logINFO, "Parameters:"); - //UBLOG(logINFO, "uLb = " << uLB); - //UBLOG(logINFO, "rho = " << rhoLB); - //UBLOG(logINFO, "nuLb = " << nuLB); - UBLOG(logINFO, "Re = " << Re); - UBLOG(logINFO, "dx = " << dx); - UBLOG(logINFO, "number of levels = " << refineLevel + 1); - UBLOG(logINFO, "numOfThreads = " << numOfThreads); - UBLOG(logINFO, "path = " << pathname); - } - - migSimulationObserver->restart((int)restartStep); - grid->setTimeStep(restartStep); - - if (myid == 0) - UBLOG(logINFO, "Restart - end"); + SPtr<UbScheduler> geoSch(new UbScheduler(1)); + SPtr<SimulationObserver> ppgeo(new WriteBoundaryConditionsSimulationObserver(grid, geoSch, pathname, WbWriterVtkXmlBinary::getInstance(), comm)); + ppgeo->update(0); + ppgeo.reset(); + } + if (myid == 0) + UBLOG(logINFO, "Preprocess - end"); + } + else + { + if (myid == 0) { + UBLOG(logINFO, "Parameters:"); + //UBLOG(logINFO, "uLb = " << uLB); + //UBLOG(logINFO, "rho = " << rhoLB); + //UBLOG(logINFO, "nuLb = " << nuLB); + UBLOG(logINFO, "Re = " << Re); + UBLOG(logINFO, "dx = " << dx); + UBLOG(logINFO, "number of levels = " << refineLevel + 1); + UBLOG(logINFO, "numOfThreads = " << numOfThreads); + UBLOG(logINFO, "path = " << pathname); } - // set connectors - SPtr<Interpolator> iProcessor(new CompressibleOffsetMomentsInterpolator()); - OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm); - grid->accept(setConnsVisitor); - - - SPtr<UbScheduler> visSch(new UbScheduler(outTime/*,beginTime,endTime*/)); - SPtr<SimulationObserver> pp(new WriteMacroscopicQuantitiesSimulationObserver(grid, visSch, pathname, WbWriterVtkXmlBinary::getInstance(), conv, comm)); - - //SPtr<UbScheduler> tavSch(new UbScheduler(100, timeAvStart, timeAvStop)); - //SPtr<TimeAveragedValuesSimulationObserver> tav(new TimeAveragedValuesSimulationObserver(grid, pathname, WbWriterVtkXmlBinary::getInstance(), tavSch, comm, - //TimeAveragedValuesSimulationObserver::Density | TimeAveragedValuesSimulationObserver::Velocity | TimeAveragedValuesSimulationObserver::Fluctuations)); - //tav->setWithGhostLayer(true); - - SPtr<UbScheduler> nupsSch(new UbScheduler(500, 1000, 3000)); - //OpenMP threads control + migSimulationObserver->restart((int)restartStep); + grid->setTimeStep(restartStep); + + if (myid == 0) + UBLOG(logINFO, "Restart - end"); + } + // set connectors + SPtr<Interpolator> iProcessor(new CompressibleOffsetMomentsInterpolator()); + OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm); + grid->accept(setConnsVisitor); + + + SPtr<UbScheduler> visSch(new UbScheduler(outTime/*,beginTime,endTime*/)); + SPtr<SimulationObserver> pp(new WriteMacroscopicQuantitiesSimulationObserver(grid, visSch, pathname, WbWriterVtkXmlBinary::getInstance(), conv, comm)); + + //SPtr<UbScheduler> tavSch(new UbScheduler(100, timeAvStart, timeAvStop)); + //SPtr<TimeAveragedValuesSimulationObserver> tav(new TimeAveragedValuesSimulationObserver(grid, pathname, WbWriterVtkXmlBinary::getInstance(), tavSch, comm, + //TimeAveragedValuesSimulationObserver::Density | TimeAveragedValuesSimulationObserver::Velocity | TimeAveragedValuesSimulationObserver::Fluctuations)); + //tav->setWithGhostLayer(true); + + + SPtr<UbScheduler> nupsSch(new UbScheduler(500, 1000, 3000)); + //OpenMP threads control #ifdef _OPENMP - omp_set_num_threads(numOfThreads); + omp_set_num_threads(numOfThreads); #endif - SPtr<SimulationObserver> npr(new NUPSCounterSimulationObserver(grid, nupsSch, numOfThreads, comm)); + SPtr<SimulationObserver> npr(new NUPSCounterSimulationObserver(grid, nupsSch, numOfThreads, comm)); - SPtr<UbScheduler> stepGhostLayer(visSch); - SPtr<Simulation> calculator(new Simulation(grid, stepGhostLayer, int(endTime))); + SPtr<UbScheduler> stepGhostLayer(visSch); + SPtr<Simulation> calculator(new Simulation(grid, stepGhostLayer, int(endTime))); - //calculator->addSimulationObserver(nupr); - calculator->addSimulationObserver(npr); - calculator->addSimulationObserver(pp); - calculator->addSimulationObserver(migSimulationObserver); - //calculator->addSimulationObserver(tav); + //calculator->addSimulationObserver(nupr); + calculator->addSimulationObserver(npr); + calculator->addSimulationObserver(pp); + calculator->addSimulationObserver(migSimulationObserver); + //calculator->addSimulationObserver(tav); - if (myid == 0) - UBLOG(logINFO, "Simulation-start"); - calculator->run(); - if (myid == 0) - UBLOG(logINFO, "Simulation-end"); - } catch (std::exception &e) { - cerr << e.what() << endl << flush; - } catch (std::string &s) { - cerr << s << endl; - } catch (...) { - cerr << "unknown exception" << endl; - } + if (myid == 0) + UBLOG(logINFO, "Simulation-start"); + calculator->run(); + if (myid == 0) + UBLOG(logINFO, "Simulation-end"); } -int main(int argc, char *argv[]) + +int main(int argc, char* argv[]) { - if (argv != NULL) { - if (argv[1] != NULL) { - run(string(argv[1])); - } else { - cout << "Configuration file is missing!" << endl; - } + try { + vf::logging::Logger::initializeLogger(); + vf::basics::ConfigurationFile config = vf::basics::loadConfig(argc, argv); + run(config); + } catch (const std::exception& e) { + VF_LOG_WARNING("{}", e.what()); + return 1; } + + return 0; } //! \} diff --git a/apps/cpu/LaminarPipeFlow/LaminarPipeFlow.cpp b/apps/cpu/LaminarPipeFlow/LaminarPipeFlow.cpp index 5bf8a37916160ec59f0db48203d82918b5a0bdc5..8277b6df195cd92fcde10189b945c27d21686bfd 100644 --- a/apps/cpu/LaminarPipeFlow/LaminarPipeFlow.cpp +++ b/apps/cpu/LaminarPipeFlow/LaminarPipeFlow.cpp @@ -39,272 +39,263 @@ using namespace std; - -void run(string configname) +void run(const vf::basics::ConfigurationFile& config) { using namespace vf::lbm::dir; - try - { - vf::basics::ConfigurationFile config; - config.load(configname); - - string pathname = config.getValue<string>("pathname"); - int numOfThreads = config.getValue<int>("numOfThreads"); - vector<int> blocknx = config.getVector<int>("blocknx"); - real endTime = config.getValue<real>("endTime"); - real outTime = config.getValue<real>("outTime"); - int refineLevel = config.getValue<int>("refineLevel"); - real dx = config.getValue<real>("dx"); - vector<real> length = config.getVector<real>("length"); - real restartStep = config.getValue<real>("restartStep"); - real cpStart = config.getValue<real>("cpStart"); - real cpStep = config.getValue<real>("cpStep"); - bool newStart = config.getValue<bool>("newStart"); - - SPtr<vf::parallel::Communicator> comm = vf::parallel::MPICommunicator::getInstance(); - int myid = comm->getProcessID(); - - real dLB = length[1] / dx; - real rhoLB1 = 0.00001; - real rhoLB2 = 0.0; - real nuLB = 0.0064; - - SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter()); - - //boundary conditions - ////////////////////////////////////////////////////////////////////////////// - SPtr<BC> noSlipBC(new NoSlipBC()); - noSlipBC->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); - - SPtr<BC> pressureBC1(new PressureBC(rhoLB1)); - pressureBC1->setBCStrategy(SPtr<BCStrategy>(new PressureNonEquilibrium())); - - SPtr<BC> pressureBC2(new PressureBC(rhoLB2)); - pressureBC2->setBCStrategy(SPtr<BCStrategy>(new PressureNonEquilibrium())); - - ////////////////////////////////////////////////////////////////////////////////// - //BC visitor - BoundaryConditionsBlockVisitor bcVisitor; - - SPtr<Grid3D> grid(new Grid3D(comm)); - - SPtr<BCSet> bcProc; - bcProc = SPtr<BCSet>(new BCSet()); - - SPtr<LBMKernel> kernel = SPtr<LBMKernel>(new K17CompressibleNavierStokes()); - kernel->setBCSet(bcProc); - - SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, d00M)); - - - ////////////////////////////////////////////////////////////////////////// - //restart - SPtr<UbScheduler> mSch(new UbScheduler(cpStep, cpStart)); - SPtr<MPIIOMigrationSimulationObserver> restart(new MPIIOMigrationSimulationObserver(grid, mSch, metisVisitor, pathname, comm)); - restart->setLBMKernel(kernel); - restart->setBCSet(bcProc); - ////////////////////////////////////////////////////////////////////////// - real R = length[1] / 2.0; - real dp = (rhoLB1 / 3. - rhoLB2 / 3.); - real L = length[0]; - real u_max = R* R / (4. * nuLB) * (dp / L); - real Re = u_max * 2 * R / nuLB; - - if (myid == 0) - { - VF_LOG_INFO("Parameters:"); - VF_LOG_INFO("p1 = {}", rhoLB1 / 3.); - VF_LOG_INFO("p1 = {}", rhoLB2 / 3.); - VF_LOG_INFO("nuLb = {}", nuLB); - VF_LOG_INFO("u_max = {}", u_max ); - VF_LOG_INFO("Re = {}", Re); - VF_LOG_INFO("dx = {}", dx); - VF_LOG_INFO("number of levels = {}", refineLevel + 1); - VF_LOG_INFO("numOfThreads = {}", numOfThreads); - VF_LOG_INFO("path = {}", pathname); - VF_LOG_INFO("Preprocess - start"); - } - - if (newStart) - { - - //bounding box - real g_minX1 = 0.0; - real g_minX2 = -length[1] / 2.0; - real g_minX3 = -length[2] / 2.0; - - real g_maxX1 = length[0]; - real g_maxX2 = length[1] / 2.0; - real g_maxX3 = length[2] / 2.0; - - SPtr<GbObject3D> cylinder(new GbCylinder3D(g_minX1 - 2.0*dx, 0.0, 0.0, g_maxX1 + 2.0*dx, 0.0, 0.0, dLB / 2.0)); - GbSystem3D::writeGeoObject(cylinder.get(), pathname + "/geo/cylinder", WbWriterVtkXmlBinary::getInstance()); - - SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3)); - if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance()); - - real blockLength = 3. * dx; - - grid->setDeltaX(dx); - grid->setBlockNX(blocknx[0], blocknx[1], blocknx[2]); - grid->setPeriodicX1(false); - grid->setPeriodicX2(false); - grid->setPeriodicX3(false); - - if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance()); - - GenBlocksGridVisitor genBlocks(gridCube); - grid->accept(genBlocks); - - SPtr<GbObject3D> refineCube1_1(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_minX2 + blockLength, g_maxX3)); - if (myid == 0) GbSystem3D::writeGeoObject(refineCube1_1.get(), pathname + "/geo/refineCube1_1", WbWriterVtkXmlBinary::getInstance()); + string pathname = config.getValue<string>("pathname"); + int numOfThreads = config.getValue<int>("numOfThreads"); + vector<int> blocknx = config.getVector<int>("blocknx"); + real endTime = config.getValue<real>("endTime"); + real outTime = config.getValue<real>("outTime"); + int refineLevel = config.getValue<int>("refineLevel"); + real dx = config.getValue<real>("dx"); + vector<real> length = config.getVector<real>("length"); + real restartStep = config.getValue<real>("restartStep"); + real cpStart = config.getValue<real>("cpStart"); + real cpStep = config.getValue<real>("cpStep"); + bool newStart = config.getValue<bool>("newStart"); + + SPtr<vf::parallel::Communicator> comm = vf::parallel::MPICommunicator::getInstance(); + int myid = comm->getProcessID(); + + real dLB = length[1] / dx; + real rhoLB1 = 0.00001; + real rhoLB2 = 0.0; + real nuLB = 0.0064; + + SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter()); + + // boundary conditions + ////////////////////////////////////////////////////////////////////////////// + SPtr<BC> noSlipBC(new NoSlipBC()); + noSlipBC->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); + + SPtr<BC> pressureBC1(new PressureBC(rhoLB1)); + pressureBC1->setBCStrategy(SPtr<BCStrategy>(new PressureNonEquilibrium())); + + SPtr<BC> pressureBC2(new PressureBC(rhoLB2)); + pressureBC2->setBCStrategy(SPtr<BCStrategy>(new PressureNonEquilibrium())); + + ////////////////////////////////////////////////////////////////////////////////// + // BC visitor + BoundaryConditionsBlockVisitor bcVisitor; + + SPtr<Grid3D> grid(new Grid3D(comm)); + + SPtr<BCSet> bcProc; + bcProc = SPtr<BCSet>(new BCSet()); + + SPtr<LBMKernel> kernel = SPtr<LBMKernel>(new K17CompressibleNavierStokes()); + kernel->setBCSet(bcProc); + + SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, d00M)); + + ////////////////////////////////////////////////////////////////////////// + // restart + SPtr<UbScheduler> mSch(new UbScheduler(cpStep, cpStart)); + SPtr<MPIIOMigrationSimulationObserver> restart( + new MPIIOMigrationSimulationObserver(grid, mSch, metisVisitor, pathname, comm)); + restart->setLBMKernel(kernel); + restart->setBCSet(bcProc); + ////////////////////////////////////////////////////////////////////////// + real R = length[1] / 2.0; + real dp = (rhoLB1 / 3. - rhoLB2 / 3.); + real L = length[0]; + real u_max = R * R / (4. * nuLB) * (dp / L); + real Re = u_max * 2 * R / nuLB; + + if (myid == 0) { + VF_LOG_INFO("Parameters:"); + VF_LOG_INFO("p1 = {}", rhoLB1 / 3.); + VF_LOG_INFO("p1 = {}", rhoLB2 / 3.); + VF_LOG_INFO("nuLb = {}", nuLB); + VF_LOG_INFO("u_max = {}", u_max); + VF_LOG_INFO("Re = {}", Re); + VF_LOG_INFO("dx = {}", dx); + VF_LOG_INFO("number of levels = {}", refineLevel + 1); + VF_LOG_INFO("numOfThreads = {}", numOfThreads); + VF_LOG_INFO("path = {}", pathname); + VF_LOG_INFO("Preprocess - start"); + } + + if (newStart) { + + // bounding box + real g_minX1 = 0.0; + real g_minX2 = -length[1] / 2.0; + real g_minX3 = -length[2] / 2.0; + + real g_maxX1 = length[0]; + real g_maxX2 = length[1] / 2.0; + real g_maxX3 = length[2] / 2.0; + + SPtr<GbObject3D> cylinder(new GbCylinder3D(g_minX1 - 2.0 * dx, 0.0, 0.0, g_maxX1 + 2.0 * dx, 0.0, 0.0, dLB / 2.0)); + GbSystem3D::writeGeoObject(cylinder.get(), pathname + "/geo/cylinder", WbWriterVtkXmlBinary::getInstance()); + + SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3)); + if (myid == 0) + GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance()); + + real blockLength = 3. * dx; + + grid->setDeltaX(dx); + grid->setBlockNX(blocknx[0], blocknx[1], blocknx[2]); + grid->setPeriodicX1(false); + grid->setPeriodicX2(false); + grid->setPeriodicX3(false); + + if (myid == 0) + GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance()); + + GenBlocksGridVisitor genBlocks(gridCube); + grid->accept(genBlocks); - SPtr<GbObject3D> refineCube1_2(new GbCuboid3D(g_minX1, g_maxX2 - blockLength, g_minX3, g_maxX1, g_maxX2, g_maxX3)); - if (myid == 0) GbSystem3D::writeGeoObject(refineCube1_2.get(), pathname + "/geo/refineCube1_2", WbWriterVtkXmlBinary::getInstance()); + SPtr<GbObject3D> refineCube1_1(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_minX2 + blockLength, g_maxX3)); + if (myid == 0) + GbSystem3D::writeGeoObject(refineCube1_1.get(), pathname + "/geo/refineCube1_1", + WbWriterVtkXmlBinary::getInstance()); - SPtr<GbObject3D> refineCube1_3(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_minX3 + blockLength)); - if (myid == 0) GbSystem3D::writeGeoObject(refineCube1_3.get(), pathname + "/geo/refineCube1_3", WbWriterVtkXmlBinary::getInstance()); + SPtr<GbObject3D> refineCube1_2(new GbCuboid3D(g_minX1, g_maxX2 - blockLength, g_minX3, g_maxX1, g_maxX2, g_maxX3)); + if (myid == 0) + GbSystem3D::writeGeoObject(refineCube1_2.get(), pathname + "/geo/refineCube1_2", + WbWriterVtkXmlBinary::getInstance()); - SPtr<GbObject3D> refineCube1_4(new GbCuboid3D(g_minX1, g_minX2, g_maxX3 - blockLength, g_maxX1, g_maxX2, g_maxX3)); - if (myid == 0) GbSystem3D::writeGeoObject(refineCube1_4.get(), pathname + "/geo/refineCube1_4", WbWriterVtkXmlBinary::getInstance()); + SPtr<GbObject3D> refineCube1_3(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_minX3 + blockLength)); + if (myid == 0) + GbSystem3D::writeGeoObject(refineCube1_3.get(), pathname + "/geo/refineCube1_3", + WbWriterVtkXmlBinary::getInstance()); - if (refineLevel > 0) - { - if (myid == 0) UBLOG(logINFO, "Refinement - start"); + SPtr<GbObject3D> refineCube1_4(new GbCuboid3D(g_minX1, g_minX2, g_maxX3 - blockLength, g_maxX1, g_maxX2, g_maxX3)); + if (myid == 0) + GbSystem3D::writeGeoObject(refineCube1_4.get(), pathname + "/geo/refineCube1_4", + WbWriterVtkXmlBinary::getInstance()); + + if (refineLevel > 0) { + if (myid == 0) + UBLOG(logINFO, "Refinement - start"); RefineCrossAndInsideGbObjectHelper refineHelper(grid, refineLevel, comm); refineHelper.addGbObject(refineCube1_1, 1); refineHelper.addGbObject(refineCube1_2, 1); refineHelper.addGbObject(refineCube1_3, 1); refineHelper.addGbObject(refineCube1_4, 1); refineHelper.refine(); - if (myid == 0) UBLOG(logINFO, "Refinement - end"); - } + if (myid == 0) + UBLOG(logINFO, "Refinement - end"); + } + + // inflow + GbCuboid3DPtr geoInflow(new GbCuboid3D(g_minX1 - blockLength, g_minX2 - blockLength, g_minX3 - blockLength, g_minX1, + g_maxX2 + blockLength, g_maxX3 + blockLength)); + if (myid == 0) + GbSystem3D::writeGeoObject(geoInflow.get(), pathname + "/geo/geoInflow", 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(), pathname + "/geo/geoOutflow", WbWriterVtkXmlASCII::getInstance()); - //inflow - GbCuboid3DPtr geoInflow(new GbCuboid3D(g_minX1-blockLength, g_minX2-blockLength, g_minX3-blockLength, g_minX1, g_maxX2+blockLength, g_maxX3+blockLength)); - if (myid==0) GbSystem3D::writeGeoObject(geoInflow.get(), pathname+"/geo/geoInflow", WbWriterVtkXmlASCII::getInstance()); + SPtr<D3Q27Interactor> cylinderInt(new D3Q27Interactor(cylinder, grid, noSlipBC, Interactor3D::INVERSESOLID)); - //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(), pathname+"/geo/geoOutflow", WbWriterVtkXmlASCII::getInstance()); + SPtr<D3Q27Interactor> inflowInt = + SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow, grid, pressureBC1, Interactor3D::SOLID)); - SPtr<D3Q27Interactor> cylinderInt(new D3Q27Interactor(cylinder, grid, noSlipBC, Interactor3D::INVERSESOLID)); - - SPtr<D3Q27Interactor> inflowInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow, grid, pressureBC1, Interactor3D::SOLID)); + SPtr<D3Q27Interactor> outflowInt = + SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow, grid, pressureBC2, Interactor3D::SOLID)); - SPtr<D3Q27Interactor> outflowInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow, grid, pressureBC2, Interactor3D::SOLID)); + InteractorsHelper intHelper(grid, metisVisitor); + intHelper.addInteractor(cylinderInt); + intHelper.addInteractor(inflowInt); + intHelper.addInteractor(outflowInt); + intHelper.selectBlocks(); - InteractorsHelper intHelper(grid, metisVisitor); - intHelper.addInteractor(cylinderInt); - intHelper.addInteractor(inflowInt); - intHelper.addInteractor(outflowInt); - intHelper.selectBlocks(); + SPtr<SimulationObserver> ppblocks(new WriteBlocksSimulationObserver( + grid, SPtr<UbScheduler>(new UbScheduler(1)), pathname, WbWriterVtkXmlBinary::getInstance(), comm)); + ppblocks->update(0); + ppblocks.reset(); - SPtr<SimulationObserver> ppblocks(new WriteBlocksSimulationObserver(grid, SPtr<UbScheduler>(new UbScheduler(1)), pathname, WbWriterVtkXmlBinary::getInstance(), comm)); - ppblocks->update(0); - ppblocks.reset(); + if (myid == 0) + VF_LOG_INFO("{}", Utilities::toString(grid, comm->getNumberOfProcesses())); - if (myid == 0) VF_LOG_INFO("{}",Utilities::toString(grid, comm->getNumberOfProcesses())); - - SetKernelBlockVisitor kernelVisitor(kernel, nuLB); - grid->accept(kernelVisitor); + SetKernelBlockVisitor kernelVisitor(kernel, nuLB); + grid->accept(kernelVisitor); - if (refineLevel > 0) - { + if (refineLevel > 0) { SetUndefinedNodesBlockVisitor undefNodesVisitor; grid->accept(undefNodesVisitor); - } + } - intHelper.setBC(); + intHelper.setBC(); - //initialization of distributions - InitDistributionsBlockVisitor initVisitor; - grid->accept(initVisitor); + // initialization of distributions + InitDistributionsBlockVisitor initVisitor; + grid->accept(initVisitor); - //boundary conditions grid - { + // boundary conditions grid + { SPtr<UbScheduler> geoSch(new UbScheduler(1)); - SPtr<SimulationObserver> ppgeo(new WriteBoundaryConditionsSimulationObserver(grid, geoSch, pathname, WbWriterVtkXmlBinary::getInstance(), comm)); + SPtr<SimulationObserver> ppgeo(new WriteBoundaryConditionsSimulationObserver( + grid, geoSch, pathname, WbWriterVtkXmlBinary::getInstance(), comm)); ppgeo->update(0); ppgeo.reset(); - } + } - - } - else - { - restart->restart((int)restartStep); - grid->setTimeStep(restartStep); + } else { + restart->restart((int)restartStep); + grid->setTimeStep(restartStep); - if (myid == 0) VF_LOG_INFO("Restart - end"); - } + if (myid == 0) + VF_LOG_INFO("Restart - end"); + } - grid->accept(bcVisitor); + grid->accept(bcVisitor); - OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm); - grid->accept(setConnsVisitor); + OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm); + grid->accept(setConnsVisitor); - SPtr<Interpolator> iProcessor(new CompressibleOffsetMomentsInterpolator()); - SetInterpolationConnectorsBlockVisitor setInterConnsVisitor(comm, nuLB, iProcessor); - grid->accept(setInterConnsVisitor); + SPtr<Interpolator> iProcessor(new CompressibleOffsetMomentsInterpolator()); + SetInterpolationConnectorsBlockVisitor setInterConnsVisitor(comm, nuLB, iProcessor); + grid->accept(setInterConnsVisitor); - SPtr<UbScheduler> visSch(new UbScheduler(outTime)); - SPtr<SimulationObserver> pp(new WriteMacroscopicQuantitiesSimulationObserver(grid, visSch, pathname, WbWriterVtkXmlBinary::getInstance(), conv, comm)); + SPtr<UbScheduler> visSch(new UbScheduler(outTime)); + SPtr<SimulationObserver> pp(new WriteMacroscopicQuantitiesSimulationObserver( + grid, visSch, pathname, WbWriterVtkXmlBinary::getInstance(), conv, comm)); - SPtr<UbScheduler> nupsSch(new UbScheduler(10, 10, 100)); - SPtr<SimulationObserver> npr(new NUPSCounterSimulationObserver(grid, nupsSch, numOfThreads, comm)); + SPtr<UbScheduler> nupsSch(new UbScheduler(10, 10, 100)); + SPtr<SimulationObserver> npr(new NUPSCounterSimulationObserver(grid, nupsSch, numOfThreads, comm)); #ifdef _OPENMP - omp_set_num_threads(numOfThreads); + omp_set_num_threads(numOfThreads); #endif - if (myid == 0) - VF_LOG_INFO("Preprocess - end"); - - SPtr<UbScheduler> stepGhostLayer(visSch); - SPtr<Simulation> simulation(new Simulation(grid, stepGhostLayer, int(endTime))); - simulation->addSimulationObserver(npr); - simulation->addSimulationObserver(pp); - simulation->addSimulationObserver(restart); - - if (myid == 0) VF_LOG_INFO("Simulation-start"); - simulation->run(); - if (myid == 0) VF_LOG_INFO("Simulation-end"); - } - catch (std::exception& e) - { - cerr << e.what() << endl << flush; - } - catch (std::string& s) - { - cerr << s << endl; - } - catch (...) - { - cerr << "unknown exception" << endl; - } + if (myid == 0) + VF_LOG_INFO("Preprocess - end"); -} + SPtr<UbScheduler> stepGhostLayer(visSch); + SPtr<Simulation> simulation(new Simulation(grid, stepGhostLayer, int(endTime))); + simulation->addSimulationObserver(npr); + simulation->addSimulationObserver(pp); + simulation->addSimulationObserver(restart); + if (myid == 0) + VF_LOG_INFO("Simulation-start"); + simulation->run(); + if (myid == 0) + VF_LOG_INFO("Simulation-end"); +} -int main(int argc, char *argv[]) +int main(int argc, char* argv[]) { try { vf::logging::Logger::initializeLogger(); - - VF_LOG_INFO("Starting VirtualFluids..."); - - if (argc > 1) - run(std::string(argv[1])); - else - VF_LOG_CRITICAL("Configuration file is missing!"); - - VF_LOG_INFO("VirtualFluids is finished."); - - } catch (const spdlog::spdlog_ex &ex) { - std::cout << "Log initialization failed: " << ex.what() << std::endl; + auto config = vf::basics::loadConfig(argc, argv); + run(config); + } catch (const std::exception& e) { + VF_LOG_WARNING("{}", e.what()); + return 1; } + return 0; } //! \} diff --git a/apps/cpu/LaminarPlaneFlow/LaminarPlaneFlow.cpp b/apps/cpu/LaminarPlaneFlow/LaminarPlaneFlow.cpp index e55145e083046e03adf5975e4e3c2c6b448263c7..7805a797ad465def68d4c5a7884914e110fb3119 100644 --- a/apps/cpu/LaminarPlaneFlow/LaminarPlaneFlow.cpp +++ b/apps/cpu/LaminarPlaneFlow/LaminarPlaneFlow.cpp @@ -39,110 +39,109 @@ using namespace std; - -void pflowdp(string configname) +void run(const vf::basics::ConfigurationFile& config) { - using namespace vf::lbm::dir; - try - { - vf::basics::ConfigurationFile config; - config.load(configname); - - string pathname = config.getValue<string>("pathname"); - int numOfThreads = config.getValue<int>("numOfThreads"); - vector<int> blocknx = config.getVector<int>("blocknx"); - vector<real> boundingBox = config.getVector<real>("boundingBox"); - real endTime = config.getValue<real>("endTime"); - real outTime = config.getValue<real>("outTime"); - int refineLevel = config.getValue<int>("refineLevel"); - real restartStep = config.getValue<real>("restartStep"); - real deltax = config.getValue<real>("deltax"); - real cpStep = config.getValue<real>("cpStep"); - real cpStart = config.getValue<real>("cpStart"); - bool newStart = config.getValue<bool>("newStart"); - - SPtr<vf::parallel::Communicator> comm = vf::parallel::MPICommunicator::getInstance(); - int myid = comm->getProcessID(); - - - real rhoLB1 = 0.00001; - real rhoLB2 = 0.0; - real nu = 0.0064; - - real h = boundingBox[1] / 2.0; - real dp = (rhoLB1 / 3. - rhoLB2 / 3.); - real L = boundingBox[0]; - real u_max = h * h / (2. * nu) * (dp / L); - real Re = u_max * 2 * h / nu; - - SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter()); - - //bounding box - real g_minX1 = 0; - real g_minX2 = 0; - real g_minX3 = 0; - - real g_maxX1 = boundingBox[0]; - real g_maxX2 = boundingBox[1]; - real g_maxX3 = boundingBox[2]; - - real blockLength = 3.0 * deltax; - - //bc - SPtr<BC> noSlipBC(new NoSlipBC()); - noSlipBC->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); - - SPtr<BC> pressureBC1(new PressureBC(rhoLB1)); - pressureBC1->setBCStrategy(SPtr<BCStrategy>(new PressureNonEquilibrium())); - - SPtr<BC> pressureBC2(new PressureBC(rhoLB2)); - pressureBC2->setBCStrategy(SPtr<BCStrategy>(new PressureNonEquilibrium())); - - //BS visitor - BoundaryConditionsBlockVisitor bcVisitor; - - SPtr<Grid3D> grid(new Grid3D(comm)); - grid->setPeriodicX1(false); - grid->setPeriodicX2(false); - grid->setPeriodicX3(true); - grid->setDeltaX(deltax); - grid->setBlockNX(blocknx[0], blocknx[1], blocknx[2]); - - SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3)); - if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance()); - - real k1 = 4; - - SPtr<GbObject3D> refineCube1_1(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2 / k1 - 1.0, g_maxX3)); - if (myid == 0) GbSystem3D::writeGeoObject(refineCube1_1.get(), pathname + "/geo/refineCube1_1", WbWriterVtkXmlBinary::getInstance()); - - SPtr<GbObject3D> refineCube1_2(new GbCuboid3D(g_minX1, g_maxX2 - g_maxX2 / k1 + 1.0, g_minX3, g_maxX1, g_maxX2, g_maxX3)); - if (myid == 0) GbSystem3D::writeGeoObject(refineCube1_2.get(), pathname + "/geo/refineCube1_2", WbWriterVtkXmlBinary::getInstance()); - - SPtr<LBMKernel> kernel; - kernel = SPtr<LBMKernel>(new K17CompressibleNavierStokes()); - - SPtr<BCSet> bcProc; - bcProc = std::make_shared<BCSet>(); - kernel->setBCSet(bcProc); - - SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, d00M)); - - ////////////////////////////////////////////////////////////////////////// - // restart - SPtr<UbScheduler> mSch(new UbScheduler(cpStep, cpStart)); - SPtr<MPIIOMigrationSimulationObserver> restart(new MPIIOMigrationSimulationObserver(grid, mSch, metisVisitor, pathname, comm)); - restart->setLBMKernel(kernel); - restart->setBCSet(bcProc); - ////////////////////////////////////////////////////////////////////////// - - if (newStart) - { - GenBlocksGridVisitor genBlocks(gridCube); - grid->accept(genBlocks); - - if (myid == 0) - { + using namespace vf::lbm::dir; + + string pathname = config.getValue<string>("pathname"); + int numOfThreads = config.getValue<int>("numOfThreads"); + vector<int> blocknx = config.getVector<int>("blocknx"); + vector<real> boundingBox = config.getVector<real>("boundingBox"); + real endTime = config.getValue<real>("endTime"); + real outTime = config.getValue<real>("outTime"); + int refineLevel = config.getValue<int>("refineLevel"); + real restartStep = config.getValue<real>("restartStep"); + real deltax = config.getValue<real>("deltax"); + real cpStep = config.getValue<real>("cpStep"); + real cpStart = config.getValue<real>("cpStart"); + bool newStart = config.getValue<bool>("newStart"); + + SPtr<vf::parallel::Communicator> comm = vf::parallel::MPICommunicator::getInstance(); + int myid = comm->getProcessID(); + + real rhoLB1 = 0.00001; + real rhoLB2 = 0.0; + real nu = 0.0064; + + real h = boundingBox[1] / 2.0; + real dp = (rhoLB1 / 3. - rhoLB2 / 3.); + real L = boundingBox[0]; + real u_max = h * h / (2. * nu) * (dp / L); + real Re = u_max * 2 * h / nu; + + SPtr<LBMUnitConverter> conv = SPtr<LBMUnitConverter>(new LBMUnitConverter()); + + // bounding box + real g_minX1 = 0; + real g_minX2 = 0; + real g_minX3 = 0; + + real g_maxX1 = boundingBox[0]; + real g_maxX2 = boundingBox[1]; + real g_maxX3 = boundingBox[2]; + + real blockLength = 3.0 * deltax; + + // bc + SPtr<BC> noSlipBC(new NoSlipBC()); + noSlipBC->setBCStrategy(SPtr<BCStrategy>(new NoSlipInterpolated())); + + SPtr<BC> pressureBC1(new PressureBC(rhoLB1)); + pressureBC1->setBCStrategy(SPtr<BCStrategy>(new PressureNonEquilibrium())); + + SPtr<BC> pressureBC2(new PressureBC(rhoLB2)); + pressureBC2->setBCStrategy(SPtr<BCStrategy>(new PressureNonEquilibrium())); + + // BS visitor + BoundaryConditionsBlockVisitor bcVisitor; + + SPtr<Grid3D> grid(new Grid3D(comm)); + grid->setPeriodicX1(false); + grid->setPeriodicX2(false); + grid->setPeriodicX3(true); + grid->setDeltaX(deltax); + grid->setBlockNX(blocknx[0], blocknx[1], blocknx[2]); + + SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3)); + if (myid == 0) + GbSystem3D::writeGeoObject(gridCube.get(), pathname + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance()); + + real k1 = 4; + + SPtr<GbObject3D> refineCube1_1(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2 / k1 - 1.0, g_maxX3)); + if (myid == 0) + GbSystem3D::writeGeoObject(refineCube1_1.get(), pathname + "/geo/refineCube1_1", + WbWriterVtkXmlBinary::getInstance()); + + SPtr<GbObject3D> refineCube1_2( + new GbCuboid3D(g_minX1, g_maxX2 - g_maxX2 / k1 + 1.0, g_minX3, g_maxX1, g_maxX2, g_maxX3)); + if (myid == 0) + GbSystem3D::writeGeoObject(refineCube1_2.get(), pathname + "/geo/refineCube1_2", + WbWriterVtkXmlBinary::getInstance()); + + SPtr<LBMKernel> kernel; + kernel = SPtr<LBMKernel>(new K17CompressibleNavierStokes()); + + SPtr<BCSet> bcProc; + bcProc = std::make_shared<BCSet>(); + kernel->setBCSet(bcProc); + + SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, d00M)); + + ////////////////////////////////////////////////////////////////////////// + // restart + SPtr<UbScheduler> mSch(new UbScheduler(cpStep, cpStart)); + SPtr<MPIIOMigrationSimulationObserver> restart( + new MPIIOMigrationSimulationObserver(grid, mSch, metisVisitor, pathname, comm)); + restart->setLBMKernel(kernel); + restart->setBCSet(bcProc); + ////////////////////////////////////////////////////////////////////////// + + if (newStart) { + GenBlocksGridVisitor genBlocks(gridCube); + grid->accept(genBlocks); + + if (myid == 0) { UBLOG(logINFO, "Parameters:"); UBLOG(logINFO, "h = " << h); UBLOG(logINFO, "nue = " << nu); @@ -154,176 +153,168 @@ void pflowdp(string configname) UBLOG(logINFO, "numOfThreads = " << numOfThreads); UBLOG(logINFO, "path = " << pathname); UBLOG(logINFO, "Preprozess - start"); - } - - //walls - GbCuboid3DPtr addWallYmin(new GbCuboid3D(g_minX1 - blockLength, g_minX2 - blockLength, g_minX3 - blockLength, g_maxX1 + blockLength, g_minX2, g_maxX3 + blockLength)); - if (myid == 0) GbSystem3D::writeGeoObject(addWallYmin.get(), pathname + "/geo/addWallYmin", WbWriterVtkXmlASCII::getInstance()); - - GbCuboid3DPtr addWallYmax(new GbCuboid3D(g_minX1 - blockLength, g_maxX2, g_minX3 - blockLength, g_maxX1 + blockLength, g_maxX2 + blockLength, g_maxX3 + blockLength)); - if (myid == 0) GbSystem3D::writeGeoObject(addWallYmax.get(), pathname + "/geo/addWallYmax", WbWriterVtkXmlASCII::getInstance()); - - //inflow - GbCuboid3DPtr geoInflow(new GbCuboid3D(g_minX1 - blockLength, g_minX2 - blockLength, g_minX3 - blockLength, g_minX1, g_maxX2 + blockLength, g_maxX3 + blockLength)); - if (myid == 0) GbSystem3D::writeGeoObject(geoInflow.get(), pathname + "/geo/geoInflow", 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(), pathname + "/geo/geoOutflow", WbWriterVtkXmlASCII::getInstance()); - - SPtr<SimulationObserver> ppblocks(new WriteBlocksSimulationObserver(grid, SPtr<UbScheduler>(new UbScheduler(1)), pathname, WbWriterVtkXmlBinary::getInstance(), comm)); - - if (refineLevel > 0) - { - if (myid == 0) UBLOG(logINFO, "Refinement - start"); + } + + // walls + GbCuboid3DPtr addWallYmin(new GbCuboid3D(g_minX1 - blockLength, g_minX2 - blockLength, g_minX3 - blockLength, + g_maxX1 + blockLength, g_minX2, g_maxX3 + blockLength)); + if (myid == 0) + GbSystem3D::writeGeoObject(addWallYmin.get(), pathname + "/geo/addWallYmin", WbWriterVtkXmlASCII::getInstance()); + + GbCuboid3DPtr addWallYmax(new GbCuboid3D(g_minX1 - blockLength, g_maxX2, g_minX3 - blockLength, + g_maxX1 + blockLength, g_maxX2 + blockLength, g_maxX3 + blockLength)); + if (myid == 0) + GbSystem3D::writeGeoObject(addWallYmax.get(), pathname + "/geo/addWallYmax", WbWriterVtkXmlASCII::getInstance()); + + // inflow + GbCuboid3DPtr geoInflow(new GbCuboid3D(g_minX1 - blockLength, g_minX2 - blockLength, g_minX3 - blockLength, g_minX1, + g_maxX2 + blockLength, g_maxX3 + blockLength)); + if (myid == 0) + GbSystem3D::writeGeoObject(geoInflow.get(), pathname + "/geo/geoInflow", 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(), pathname + "/geo/geoOutflow", WbWriterVtkXmlASCII::getInstance()); + + SPtr<SimulationObserver> ppblocks(new WriteBlocksSimulationObserver( + grid, SPtr<UbScheduler>(new UbScheduler(1)), pathname, WbWriterVtkXmlBinary::getInstance(), comm)); + + if (refineLevel > 0) { + if (myid == 0) + UBLOG(logINFO, "Refinement - start"); RefineCrossAndInsideGbObjectHelper refineHelper(grid, refineLevel, comm); refineHelper.addGbObject(refineCube1_1, 1); refineHelper.addGbObject(refineCube1_2, 1); refineHelper.refine(); - if (myid == 0) UBLOG(logINFO, "Refinement - end"); - } - - //walls - SPtr<D3Q27Interactor> addWallYminInt(new D3Q27Interactor(addWallYmin, grid, noSlipBC, Interactor3D::SOLID)); - SPtr<D3Q27Interactor> addWallYmaxInt(new D3Q27Interactor(addWallYmax, grid, noSlipBC, Interactor3D::SOLID)); + if (myid == 0) + UBLOG(logINFO, "Refinement - end"); + } - //inflow - SPtr<D3Q27Interactor> inflowInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow, grid, pressureBC1, Interactor3D::SOLID)); + // walls + SPtr<D3Q27Interactor> addWallYminInt(new D3Q27Interactor(addWallYmin, grid, noSlipBC, Interactor3D::SOLID)); + SPtr<D3Q27Interactor> addWallYmaxInt(new D3Q27Interactor(addWallYmax, grid, noSlipBC, Interactor3D::SOLID)); - //outflow - SPtr<D3Q27Interactor> outflowInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow, grid, pressureBC2, Interactor3D::SOLID)); - //SPtr<D3Q27Interactor> outflowSolidInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow, grid, noSlipBC, Interactor3D::SOLID)); + // inflow + SPtr<D3Q27Interactor> inflowInt = + SPtr<D3Q27Interactor>(new D3Q27Interactor(geoInflow, grid, pressureBC1, Interactor3D::SOLID)); - //////////////////////////////////////////// - //METIS - SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, d00M)); - //////////////////////////////////////////// - /////delete solid blocks - if (myid == 0) UBLOG(logINFO, "deleteSolidBlocks - start"); - InteractorsHelper intHelper(grid, metisVisitor); + // outflow + SPtr<D3Q27Interactor> outflowInt = + SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow, grid, pressureBC2, Interactor3D::SOLID)); + // SPtr<D3Q27Interactor> outflowSolidInt = SPtr<D3Q27Interactor>(new D3Q27Interactor(geoOutflow, grid, noSlipBC, + // Interactor3D::SOLID)); - intHelper.addInteractor(addWallYminInt); - intHelper.addInteractor(addWallYmaxInt); + //////////////////////////////////////////// + // METIS + SPtr<Grid3DVisitor> metisVisitor( + new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, d00M)); + //////////////////////////////////////////// + /////delete solid blocks + if (myid == 0) + UBLOG(logINFO, "deleteSolidBlocks - start"); + InteractorsHelper intHelper(grid, metisVisitor); - intHelper.addInteractor(inflowInt); + intHelper.addInteractor(addWallYminInt); + intHelper.addInteractor(addWallYmaxInt); - intHelper.addInteractor(outflowInt); + intHelper.addInteractor(inflowInt); + intHelper.addInteractor(outflowInt); + intHelper.selectBlocks(); + if (myid == 0) + UBLOG(logINFO, "deleteSolidBlocks - end"); + ////////////////////////////////////// - intHelper.selectBlocks(); - if (myid == 0) UBLOG(logINFO, "deleteSolidBlocks - end"); - ////////////////////////////////////// + ppblocks->update(0); + ppblocks.reset(); + if (myid == 0) + VF_LOG_INFO("{}", Utilities::toString(grid, comm->getNumberOfProcesses())); + SetKernelBlockVisitor kernelVisitor(kernel, nu); + grid->accept(kernelVisitor); - ppblocks->update(0); - ppblocks.reset(); - - if (myid == 0) VF_LOG_INFO("{}",Utilities::toString(grid, comm->getNumberOfProcesses())); - - - SetKernelBlockVisitor kernelVisitor(kernel, nu); - grid->accept(kernelVisitor); - - if (refineLevel > 0) - { + if (refineLevel > 0) { SetUndefinedNodesBlockVisitor undefNodesVisitor; grid->accept(undefNodesVisitor); - } - - //walls - intHelper.setBC(); + } - grid->accept(bcVisitor); + // walls + intHelper.setBC(); + grid->accept(bcVisitor); - InitDistributionsBlockVisitor initVisitor; - grid->accept(initVisitor); + InitDistributionsBlockVisitor initVisitor; + grid->accept(initVisitor); - //Postprocess - SPtr<UbScheduler> geoSch(new UbScheduler(1)); - SPtr<SimulationObserver> ppgeo( - new WriteBoundaryConditionsSimulationObserver(grid, geoSch, pathname, WbWriterVtkXmlBinary::getInstance(), comm)); - ppgeo->update(0); - ppgeo.reset(); + // Postprocess + SPtr<UbScheduler> geoSch(new UbScheduler(1)); + SPtr<SimulationObserver> ppgeo(new WriteBoundaryConditionsSimulationObserver( + grid, geoSch, pathname, WbWriterVtkXmlBinary::getInstance(), comm)); + ppgeo->update(0); + ppgeo.reset(); - if (myid == 0) UBLOG(logINFO, "Preprocess - end"); - } - else - { - restart->restart((int)restartStep); - grid->setTimeStep(restartStep); + if (myid == 0) + UBLOG(logINFO, "Preprocess - end"); + } else { + restart->restart((int)restartStep); + grid->setTimeStep(restartStep); - if (myid == 0) UBLOG(logINFO, "Restart - end"); - } + if (myid == 0) + UBLOG(logINFO, "Restart - end"); + } - grid->accept(bcVisitor); + grid->accept(bcVisitor); - // set connectors - OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm); - grid->accept(setConnsVisitor); + // set connectors + OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm); + grid->accept(setConnsVisitor); - SPtr<Interpolator> iProcessor(new CompressibleOffsetMomentsInterpolator()); - SetInterpolationConnectorsBlockVisitor setInterConnsVisitor(comm, nu, iProcessor); - grid->accept(setInterConnsVisitor); + SPtr<Interpolator> iProcessor(new CompressibleOffsetMomentsInterpolator()); + SetInterpolationConnectorsBlockVisitor setInterConnsVisitor(comm, nu, iProcessor); + grid->accept(setInterConnsVisitor); - SPtr<UbScheduler> nupsSch(new UbScheduler(10, 30, 100)); - SPtr<SimulationObserver> npr(new NUPSCounterSimulationObserver (grid, nupsSch, numOfThreads, comm)); + SPtr<UbScheduler> nupsSch(new UbScheduler(10, 30, 100)); + SPtr<SimulationObserver> npr(new NUPSCounterSimulationObserver(grid, nupsSch, numOfThreads, comm)); - //write data for visualization of macroscopic quantities - SPtr<UbScheduler> visSch(new UbScheduler(outTime)); - SPtr<WriteMacroscopicQuantitiesSimulationObserver> writeMQSimulationObserver(new WriteMacroscopicQuantitiesSimulationObserver(grid, visSch, pathname, - WbWriterVtkXmlASCII::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm)); + // write data for visualization of macroscopic quantities + SPtr<UbScheduler> visSch(new UbScheduler(outTime)); + SPtr<WriteMacroscopicQuantitiesSimulationObserver> writeMQSimulationObserver( + new WriteMacroscopicQuantitiesSimulationObserver(grid, visSch, pathname, WbWriterVtkXmlASCII::getInstance(), + SPtr<LBMUnitConverter>(new LBMUnitConverter()), comm)); - - // OpenMP threads control + // OpenMP threads control #ifdef _OPENMP - omp_set_num_threads(numOfThreads); + omp_set_num_threads(numOfThreads); #endif - // Create simulation - SPtr<Simulation> simulation(new Simulation(grid, visSch, endTime)); - simulation->addSimulationObserver(npr); - simulation->addSimulationObserver(writeMQSimulationObserver); - - ////////////////////////////////////////////////////////////////////////// - // Run simulation - ////////////////////////////////////////////////////////////////////////// - - UBLOG(logINFO, "Simulation-start"); - simulation->run(); - UBLOG(logINFO, "Simulation-end"); - } - catch (std::exception & e) - { - cerr << e.what() << endl << flush; - } - catch (std::string & s) - { - cerr << s << endl; - } - catch (...) - { - cerr << "unknown exception" << endl; - } - + // Create simulation + SPtr<Simulation> simulation(new Simulation(grid, visSch, endTime)); + simulation->addSimulationObserver(npr); + simulation->addSimulationObserver(writeMQSimulationObserver); + + ////////////////////////////////////////////////////////////////////////// + // Run simulation + ////////////////////////////////////////////////////////////////////////// + + UBLOG(logINFO, "Simulation-start"); + simulation->run(); + UBLOG(logINFO, "Simulation-end"); } -////////////////////////////////////////////////////////////////////////// + int main(int argc, char* argv[]) { - if (argv != NULL) - { - if (argv[1] != NULL) - { - pflowdp(string(argv[1])); - } - else - { - cout << "Configuration file is missing!" << endl; - } - } - - return 0; + try { + vf::logging::Logger::initializeLogger(); + auto config = vf::basics::loadConfig(argc, argv); + run(config); + } catch (const std::exception& e) { + VF_LOG_WARNING("{}", e.what()); + return 1; + } + return 0; } //! \} diff --git a/apps/cpu/LidDrivenCavity/LidDrivenCavity.cpp b/apps/cpu/LidDrivenCavity/LidDrivenCavity.cpp index 98be72b2084549cd6d81e105273444ac9050f100..37c29df6de2d0c4b9d418cf114e2f94fbc26b6bf 100644 --- a/apps/cpu/LidDrivenCavity/LidDrivenCavity.cpp +++ b/apps/cpu/LidDrivenCavity/LidDrivenCavity.cpp @@ -38,228 +38,209 @@ using namespace std; -int main(int argc, char* argv[]) +int main(int argc, char* argv[]) { - using namespace vf::lbm::dir; - string configname; - if (argv != NULL) - { - if (argv[1] != NULL) - { - configname = string(argv[1]); - } - else - { - cout << "Configuration file is missing!" << endl; - return 0; - } - } - - try - { - ////////////////////////////////////////////////////////////////////////// - // Simulation parameters - ////////////////////////////////////////////////////////////////////////// - vf::basics::ConfigurationFile config; - config.load(configname); - - // set your output path here - string path = config.getValue<string>("path"); - - const real L = config.getValue<real>("L"); - const real Re = config.getValue<real>("Re"); - const real velocity = config.getValue<real>("velocity"); - const real dt = config.getValue<real>("dt"); - const unsigned int nx = config.getValue<unsigned int>("nx"); - - const real timeStepOut = config.getValue<real>("timeStepOut"); - const real timeStepEnd = config.getValue<real>("timeStepEnd"); - - // Number of OpenMP threads - int numOfThreads = config.getValue<int>("numOfThreads"); - - ////////////////////////////////////////////////////////////////////////// - - real dx = L / real(nx); - const real velocityLB = velocity * dt / dx; // LB units - const real u = velocityLB / sqrt(2.0); // LB units - const real viscosityLB = nx * velocityLB / Re; // LB unit - - ////////////////////////////////////////////////////////////////////////// - // create grid - ////////////////////////////////////////////////////////////////////////// - // bounding box - real g_minX1 = -0.5; - real g_minX2 = -0.5; - real g_minX3 = -0.5; - - real g_maxX1 = 0.5; - real g_maxX2 = 0.5; - real g_maxX3 = 0.5; - - SPtr<vf::parallel::Communicator> comm = vf::parallel::MPICommunicator::getInstance(); - int myid = comm->getProcessID(); - - // new grid object - SPtr<Grid3D> grid(new Grid3D(comm)); - // set grid spacing - grid->setDeltaX(dx); - // set block size for three dimensions - int blockSize = nx / 2; - grid->setBlockNX(blockSize,blockSize,blockSize); - - // Create simulation bounding box - SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3)); - GbSystem3D::writeGeoObject(gridCube.get(), path + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance()); - - UBLOG(logINFO, "Lid Driven Cavity:"); - UBLOG(logINFO, "Domain size = " << nx << " x "<< nx << " x "<< nx); - UBLOG(logINFO, "Block size = " << blockSize << " x "<< blockSize << " x "<< blockSize); - UBLOG(logINFO, "velocity = " << velocity << " m/s"); - UBLOG(logINFO, "velocityLB = " << velocityLB); - UBLOG(logINFO, "viscosityLB = " << viscosityLB); - UBLOG(logINFO, "u = " << u); - UBLOG(logINFO, "Re = " << Re); - UBLOG(logINFO, "dx = " << dx); - UBLOG(logINFO, "dt = " << dt); - UBLOG(logINFO, "Preprocess - start"); - - // Generate block grid - GenBlocksGridVisitor genBlocks(gridCube); - grid->accept(genBlocks); - - // Write block grid to VTK-file - auto ppblocks = std::make_shared<WriteBlocksSimulationObserver>(grid, SPtr<UbScheduler>(new UbScheduler(1)), path, WbWriterVtkXmlBinary::getInstance(), comm); - ppblocks->update(0); - ppblocks.reset(); - - // Create LBM kernel - auto kernel = std::make_shared<K17CompressibleNavierStokes>(); - - ////////////////////////////////////////////////////////////////////////// - // Create boundary conditions (BC) - ////////////////////////////////////////////////////////////////////////// - // Create no-slip BC - auto noSlipBC = std::make_shared<NoSlipBC>(); - noSlipBC->setBCStrategy(std::make_shared<NoSlipInterpolated>()); - - // Velocity BC - mu::Parser fct; - fct.SetExpr("u"); - fct.DefineConst("u", u); - // Set the same velocity in x and y-direction - auto velBC = std::make_shared<VelocityBC>(true, true, false, fct, 0, BCFunction::INFCONST); - velBC->setBCStrategy(std::make_shared<VelocityInterpolated>()); - - // Add velocity boundary condition to visitor. No-slip boundary - BoundaryConditionsBlockVisitor bcVisitor; - - // Create boundary conditions - SPtr<BCSet> bcProc; - bcProc = std::make_shared<BCSet>(); - kernel->setBCSet(bcProc); - - // Create boundary conditions geometry - GbCuboid3DPtr wallXmin(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_minX1, g_maxX2 + dx, g_maxX3)); - GbSystem3D::writeGeoObject(wallXmin.get(), path + "/geo/wallXmin", WbWriterVtkXmlASCII::getInstance()); - GbCuboid3DPtr wallXmax(new GbCuboid3D(g_maxX1, g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_maxX3)); - GbSystem3D::writeGeoObject(wallXmax.get(), path + "/geo/wallXmax", WbWriterVtkXmlASCII::getInstance()); - GbCuboid3DPtr wallYmin(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_minX2, g_maxX3)); - GbSystem3D::writeGeoObject(wallYmin.get(), path + "/geo/wallYmin", WbWriterVtkXmlASCII::getInstance()); - GbCuboid3DPtr wallYmax(new GbCuboid3D(g_minX1 - dx, g_maxX2, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_maxX3)); - GbSystem3D::writeGeoObject(wallYmax.get(), path + "/geo/wallYmax", WbWriterVtkXmlASCII::getInstance()); - GbCuboid3DPtr wallZmin(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_minX3)); - GbSystem3D::writeGeoObject(wallZmin.get(), path + "/geo/wallZmin", WbWriterVtkXmlASCII::getInstance()); - GbCuboid3DPtr wallZmax(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_maxX3, g_maxX1 + dx, g_maxX2 + dx, g_maxX3 + dx)); - GbSystem3D::writeGeoObject(wallZmax.get(), path + "/geo/wallZmax", WbWriterVtkXmlASCII::getInstance()); - - // Add boundary conditions to grid generator - SPtr<D3Q27Interactor> wallXminInt(new D3Q27Interactor(wallXmin, grid, noSlipBC, Interactor3D::SOLID)); - SPtr<D3Q27Interactor> wallXmaxInt(new D3Q27Interactor(wallXmax, grid, noSlipBC, Interactor3D::SOLID)); - SPtr<D3Q27Interactor> wallYminInt(new D3Q27Interactor(wallYmin, grid, noSlipBC, Interactor3D::SOLID)); - SPtr<D3Q27Interactor> wallYmaxInt(new D3Q27Interactor(wallYmax, grid, noSlipBC, Interactor3D::SOLID)); - SPtr<D3Q27Interactor> wallZminInt(new D3Q27Interactor(wallZmin, grid, noSlipBC, Interactor3D::SOLID)); - SPtr<D3Q27Interactor> wallZmaxInt(new D3Q27Interactor(wallZmax, grid, velBC, Interactor3D::SOLID)); - - SPtr<Grid3DVisitor> metisVisitor(new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, d00M)); - InteractorsHelper intHelper(grid, metisVisitor); - intHelper.addInteractor(wallZmaxInt); - intHelper.addInteractor(wallXminInt); - intHelper.addInteractor(wallXmaxInt); - intHelper.addInteractor(wallZminInt); - intHelper.addInteractor(wallYminInt); - intHelper.addInteractor(wallYmaxInt); - - intHelper.selectBlocks(); - - if (myid == 0) VF_LOG_INFO("{}",Utilities::toString(grid, comm->getNumberOfProcesses())); - - // Generate grid - SetKernelBlockVisitor kernelVisitor(kernel, viscosityLB); - grid->accept(kernelVisitor); - - intHelper.setBC(); - - // Initialization of distributions - InitDistributionsBlockVisitor initVisitor; - grid->accept(initVisitor); - - // Set connectivity between blocks - OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm); - grid->accept(setConnsVisitor); - - // Create lists of boundary nodes - grid->accept(bcVisitor); - - // Write grid with boundary conditions information to VTK-file - SPtr<UbScheduler> geoSch(new UbScheduler(1)); - WriteBoundaryConditionsSimulationObserver ppgeo(grid, geoSch, path, WbWriterVtkXmlBinary::getInstance(), comm); - ppgeo.update(0); - - UBLOG(logINFO, "Preprocess - end"); - - UBLOG(logINFO, "Total Physical Memory (RAM): " << Utilities::getTotalPhysMem()/1e9 << " GB"); - UBLOG(logINFO, "Physical Memory currently used: " << Utilities::getPhysMemUsed()/1e9 << " GB"); - UBLOG(logINFO, "Physical Memory currently used by current process: " << Utilities::getPhysMemUsedByMe()/1e9 << " GB"); - - // Create coprocessor object for writing macroscopic quantities to VTK-file - SPtr<UbScheduler> visSch(new UbScheduler(timeStepOut)); - SPtr<SimulationObserver> mqSimulationObserver(new WriteMacroscopicQuantitiesSimulationObserver(grid, visSch, path, WbWriterVtkXmlBinary::getInstance(), SPtr<LBMUnitConverter>(new LBMUnitConverter(L, velocity, 1.0, nx, velocityLB)), comm)); - mqSimulationObserver->update(0); - - // Create coprocessor object for writing NUPS - SPtr<UbScheduler> nupsSch(new UbScheduler(100, 100)); - SPtr<SimulationObserver> nupsSimulationObserver(new NUPSCounterSimulationObserver(grid, nupsSch, numOfThreads, comm)); - - // OpenMP threads control + using namespace vf::lbm::dir; + try { + vf::logging::Logger::initializeLogger(); + ////////////////////////////////////////////////////////////////////////// + // Simulation parameters + ////////////////////////////////////////////////////////////////////////// + vf::basics::ConfigurationFile config = vf::basics::loadConfig(argc, argv); + + // set your output path here + string path = config.getValue<string>("path"); + const real L = config.getValue<real>("L"); + const real Re = config.getValue<real>("Re"); + const real velocity = config.getValue<real>("velocity"); + const real dt = config.getValue<real>("dt"); + const unsigned int nx = config.getValue<unsigned int>("nx"); + + const real timeStepOut = config.getValue<real>("timeStepOut"); + const real timeStepEnd = config.getValue<real>("timeStepEnd"); + + // Number of OpenMP threads + int numOfThreads = config.getValue<int>("numOfThreads"); + ////////////////////////////////////////////////////////////////////////// + real dx = L / real(nx); + const real velocityLB = velocity * dt / dx; // LB units + const real u = velocityLB / sqrt(2.0); // LB units + const real viscosityLB = nx * velocityLB / Re; // LB unit + ////////////////////////////////////////////////////////////////////////// + // create grid + ////////////////////////////////////////////////////////////////////////// + // bounding box + real g_minX1 = -0.5; + real g_minX2 = -0.5; + real g_minX3 = -0.5; + + real g_maxX1 = 0.5; + real g_maxX2 = 0.5; + real g_maxX3 = 0.5; + + SPtr<vf::parallel::Communicator> comm = vf::parallel::MPICommunicator::getInstance(); + int myid = comm->getProcessID(); + + // new grid object + SPtr<Grid3D> grid(new Grid3D(comm)); + // set grid spacing + grid->setDeltaX(dx); + // set block size for three dimensions + int blockSize = nx / 2; + grid->setBlockNX(blockSize, blockSize, blockSize); + + // Create simulation bounding box + SPtr<GbObject3D> gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3)); + GbSystem3D::writeGeoObject(gridCube.get(), path + "/geo/gridCube", WbWriterVtkXmlBinary::getInstance()); + + UBLOG(logINFO, "Lid Driven Cavity:"); + UBLOG(logINFO, "Domain size = " << nx << " x " << nx << " x " << nx); + UBLOG(logINFO, "Block size = " << blockSize << " x " << blockSize << " x " << blockSize); + UBLOG(logINFO, "velocity = " << velocity << " m/s"); + UBLOG(logINFO, "velocityLB = " << velocityLB); + UBLOG(logINFO, "viscosityLB = " << viscosityLB); + UBLOG(logINFO, "u = " << u); + UBLOG(logINFO, "Re = " << Re); + UBLOG(logINFO, "dx = " << dx); + UBLOG(logINFO, "dt = " << dt); + UBLOG(logINFO, "Preprocess - start"); + + // Generate block grid + GenBlocksGridVisitor genBlocks(gridCube); + grid->accept(genBlocks); + + // Write block grid to VTK-file + auto ppblocks = std::make_shared<WriteBlocksSimulationObserver>(grid, SPtr<UbScheduler>(new UbScheduler(1)), path, + WbWriterVtkXmlBinary::getInstance(), comm); + ppblocks->update(0); + ppblocks.reset(); + + // Create LBM kernel + auto kernel = std::make_shared<K17CompressibleNavierStokes>(); + + ////////////////////////////////////////////////////////////////////////// + // Create boundary conditions (BC) + ////////////////////////////////////////////////////////////////////////// + // Create no-slip BC + auto noSlipBC = std::make_shared<NoSlipBC>(); + noSlipBC->setBCStrategy(std::make_shared<NoSlipInterpolated>()); + + // Velocity BC + mu::Parser fct; + fct.SetExpr("u"); + fct.DefineConst("u", u); + // Set the same velocity in x and y-direction + auto velBC = std::make_shared<VelocityBC>(true, true, false, fct, 0, BCFunction::INFCONST); + velBC->setBCStrategy(std::make_shared<VelocityInterpolated>()); + + // Add velocity boundary condition to visitor. No-slip boundary + BoundaryConditionsBlockVisitor bcVisitor; + + // Create boundary conditions + SPtr<BCSet> bcProc; + bcProc = std::make_shared<BCSet>(); + kernel->setBCSet(bcProc); + + // Create boundary conditions geometry + GbCuboid3DPtr wallXmin(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_minX1, g_maxX2 + dx, g_maxX3)); + GbSystem3D::writeGeoObject(wallXmin.get(), path + "/geo/wallXmin", WbWriterVtkXmlASCII::getInstance()); + GbCuboid3DPtr wallXmax(new GbCuboid3D(g_maxX1, g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_maxX3)); + GbSystem3D::writeGeoObject(wallXmax.get(), path + "/geo/wallXmax", WbWriterVtkXmlASCII::getInstance()); + GbCuboid3DPtr wallYmin(new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_minX2, g_maxX3)); + GbSystem3D::writeGeoObject(wallYmin.get(), path + "/geo/wallYmin", WbWriterVtkXmlASCII::getInstance()); + GbCuboid3DPtr wallYmax(new GbCuboid3D(g_minX1 - dx, g_maxX2, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_maxX3)); + GbSystem3D::writeGeoObject(wallYmax.get(), path + "/geo/wallYmax", WbWriterVtkXmlASCII::getInstance()); + GbCuboid3DPtr wallZmin( + new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_minX3 - dx, g_maxX1 + dx, g_maxX2 + dx, g_minX3)); + GbSystem3D::writeGeoObject(wallZmin.get(), path + "/geo/wallZmin", WbWriterVtkXmlASCII::getInstance()); + GbCuboid3DPtr wallZmax( + new GbCuboid3D(g_minX1 - dx, g_minX2 - dx, g_maxX3, g_maxX1 + dx, g_maxX2 + dx, g_maxX3 + dx)); + GbSystem3D::writeGeoObject(wallZmax.get(), path + "/geo/wallZmax", WbWriterVtkXmlASCII::getInstance()); + + // Add boundary conditions to grid generator + SPtr<D3Q27Interactor> wallXminInt(new D3Q27Interactor(wallXmin, grid, noSlipBC, Interactor3D::SOLID)); + SPtr<D3Q27Interactor> wallXmaxInt(new D3Q27Interactor(wallXmax, grid, noSlipBC, Interactor3D::SOLID)); + SPtr<D3Q27Interactor> wallYminInt(new D3Q27Interactor(wallYmin, grid, noSlipBC, Interactor3D::SOLID)); + SPtr<D3Q27Interactor> wallYmaxInt(new D3Q27Interactor(wallYmax, grid, noSlipBC, Interactor3D::SOLID)); + SPtr<D3Q27Interactor> wallZminInt(new D3Q27Interactor(wallZmin, grid, noSlipBC, Interactor3D::SOLID)); + SPtr<D3Q27Interactor> wallZmaxInt(new D3Q27Interactor(wallZmax, grid, velBC, Interactor3D::SOLID)); + + SPtr<Grid3DVisitor> metisVisitor( + new MetisPartitioningGridVisitor(comm, MetisPartitioningGridVisitor::LevelBased, d00M)); + InteractorsHelper intHelper(grid, metisVisitor); + intHelper.addInteractor(wallZmaxInt); + intHelper.addInteractor(wallXminInt); + intHelper.addInteractor(wallXmaxInt); + intHelper.addInteractor(wallZminInt); + intHelper.addInteractor(wallYminInt); + intHelper.addInteractor(wallYmaxInt); + + intHelper.selectBlocks(); + + if (myid == 0) + VF_LOG_INFO("{}", Utilities::toString(grid, comm->getNumberOfProcesses())); + + // Generate grid + SetKernelBlockVisitor kernelVisitor(kernel, viscosityLB); + grid->accept(kernelVisitor); + + intHelper.setBC(); + + // Initialization of distributions + InitDistributionsBlockVisitor initVisitor; + grid->accept(initVisitor); + + // Set connectivity between blocks + OneDistributionSetConnectorsBlockVisitor setConnsVisitor(comm); + grid->accept(setConnsVisitor); + + // Create lists of boundary nodes + grid->accept(bcVisitor); + + // Write grid with boundary conditions information to VTK-file + SPtr<UbScheduler> geoSch(new UbScheduler(1)); + WriteBoundaryConditionsSimulationObserver ppgeo(grid, geoSch, path, WbWriterVtkXmlBinary::getInstance(), comm); + ppgeo.update(0); + + UBLOG(logINFO, "Preprocess - end"); + + UBLOG(logINFO, "Total Physical Memory (RAM): " << Utilities::getTotalPhysMem() / 1e9 << " GB"); + UBLOG(logINFO, "Physical Memory currently used: " << Utilities::getPhysMemUsed() / 1e9 << " GB"); + UBLOG(logINFO, + "Physical Memory currently used by current process: " << Utilities::getPhysMemUsedByMe() / 1e9 << " GB"); + + // Create coprocessor object for writing macroscopic quantities to VTK-file + SPtr<UbScheduler> visSch(new UbScheduler(timeStepOut)); + SPtr<SimulationObserver> mqSimulationObserver(new WriteMacroscopicQuantitiesSimulationObserver( + grid, visSch, path, WbWriterVtkXmlBinary::getInstance(), + SPtr<LBMUnitConverter>(new LBMUnitConverter(L, velocity, 1.0, nx, velocityLB)), comm)); + mqSimulationObserver->update(0); + + // Create coprocessor object for writing NUPS + SPtr<UbScheduler> nupsSch(new UbScheduler(100, 100)); + SPtr<SimulationObserver> nupsSimulationObserver( + new NUPSCounterSimulationObserver(grid, nupsSch, numOfThreads, comm)); + + // OpenMP threads control #ifdef _OPENMP - omp_set_num_threads(numOfThreads); + omp_set_num_threads(numOfThreads); #endif - // Create simulation - SPtr<Simulation> simulation(new Simulation(grid, visSch, timeStepEnd)); - simulation->addSimulationObserver(nupsSimulationObserver); - simulation->addSimulationObserver(mqSimulationObserver); - - ////////////////////////////////////////////////////////////////////////// - // Run simulation - ////////////////////////////////////////////////////////////////////////// - - UBLOG(logINFO, "Simulation-start"); - simulation->run(); - UBLOG(logINFO, "Simulation-end"); - } - catch (std::exception& e) - { - cerr << e.what() << endl << flush; - } - catch (std::string& s) - { - cerr << s << endl; - } - catch (...) - { - cerr << "unknown exception" << endl; - } + // Create simulation + SPtr<Simulation> simulation(new Simulation(grid, visSch, timeStepEnd)); + simulation->addSimulationObserver(nupsSimulationObserver); + simulation->addSimulationObserver(mqSimulationObserver); + + ////////////////////////////////////////////////////////////////////////// + // Run simulation + ////////////////////////////////////////////////////////////////////////// + UBLOG(logINFO, "Simulation-start"); + simulation->run(); + UBLOG(logINFO, "Simulation-end"); + } catch (const std::exception& e) { + VF_LOG_WARNING("{}", e.what()); + return 1; + } + return 0; } //! \} diff --git a/apps/gpu/ActuatorLine/ActuatorLine.cpp b/apps/gpu/ActuatorLine/ActuatorLine.cpp index 2aebd6d632fcecfae32d89d51b25b201c8d45413..f54f34175573bf2f4c772883d6feace6780079e6 100644 --- a/apps/gpu/ActuatorLine/ActuatorLine.cpp +++ b/apps/gpu/ActuatorLine/ActuatorLine.cpp @@ -73,10 +73,8 @@ ////////////////////////////////////////////////////////////////////////// -void run(vf::basics::ConfigurationFile& config) +void run(const vf::basics::ConfigurationFile& config) { - vf::logging::Logger::initializeLogger(); - ////////////////////////////////////////////////////////////////////////// // Simulation parameters ////////////////////////////////////////////////////////////////////////// @@ -306,20 +304,13 @@ void run(vf::basics::ConfigurationFile& config) int main(int argc, char* argv[]) { - if (argv == NULL) - return 0; - try { + vf::logging::Logger::initializeLogger(); auto config = vf::basics::loadConfig(argc, argv, "./actuatorline.cfg"); run(config); - } catch (const spdlog::spdlog_ex& ex) { - std::cout << "Log initialization failed: " << ex.what() << std::endl; - } catch (const std::bad_alloc& e) { - VF_LOG_CRITICAL("Bad Alloc: {}", e.what()); } catch (const std::exception& e) { - VF_LOG_CRITICAL("exception: {}", e.what()); - } catch (...) { - VF_LOG_CRITICAL("Unknown exception!"); + VF_LOG_WARNING("{}", e.what()); + return 1; } return 0; } diff --git a/apps/gpu/DrivenCavity/DrivenCavity.cpp b/apps/gpu/DrivenCavity/DrivenCavity.cpp index cf5e194c231c1354b875b2babcbb091656a851fa..1fe5b287175663a0417abbbdafe7629741f42ef7 100644 --- a/apps/gpu/DrivenCavity/DrivenCavity.cpp +++ b/apps/gpu/DrivenCavity/DrivenCavity.cpp @@ -204,16 +204,10 @@ int main(int argc, char* argv[]) Simulation sim(para, cudaMemoryManager, communicator, *gridGenerator, &bcFactory, &scalingFactory); sim.run(); - } catch (const spdlog::spdlog_ex& ex) { - std::cout << "Log initialization failed: " << ex.what() << std::endl; - } catch (const std::bad_alloc& e) { - VF_LOG_CRITICAL("Bad Alloc: {}", e.what()); } catch (const std::exception& e) { - VF_LOG_CRITICAL("exception: {}", e.what()); - } catch (...) { - VF_LOG_CRITICAL("Unknown exception!"); + VF_LOG_WARNING("{}", e.what()); + return 1; } - return 0; } diff --git a/apps/gpu/DrivenCavityMultiGPU/DrivenCavityMultiGPU.cpp b/apps/gpu/DrivenCavityMultiGPU/DrivenCavityMultiGPU.cpp index 1e425425857ca420e7a8a61cf121851b4bec8b9d..2b91fe4aefc0ca007c3e0ce4740a29b2cc63c298 100755 --- a/apps/gpu/DrivenCavityMultiGPU/DrivenCavityMultiGPU.cpp +++ b/apps/gpu/DrivenCavityMultiGPU/DrivenCavityMultiGPU.cpp @@ -201,36 +201,22 @@ void runVirtualFluids(const vf::basics::ConfigurationFile &config) auto cudaMemoryManager = std::make_shared<CudaMemoryManager>(para); SPtr<GridProvider> gridGenerator = GridProvider::makeGridGenerator(gridBuilderFacade->getGridBuilder(), para, cudaMemoryManager, communicator); Simulation sim(para, cudaMemoryManager, communicator, *gridGenerator, &bcFactory, &scalingFactory); - - // run simulation sim.run(); } -int main(int argc, char *argv[]) +int main(int argc, char* argv[]) { - MPI_Init(&argc, &argv); - std::string str, str2, configFile; - - if (argv != NULL) { - - try { - VF_LOG_TRACE("For the default config path to work, execute the app from the project root."); - vf::basics::ConfigurationFile config = vf::basics::loadConfig(argc, argv, "./apps/gpu/DrivenCavityMultiGPU/drivencavity_1gpu.cfg"); - runVirtualFluids(config); - - ////////////////////////////////////////////////////////////////////////// - } catch (const spdlog::spdlog_ex &ex) { - std::cout << "Log initialization failed: " << ex.what() << std::endl; - } catch (const std::bad_alloc &e) { - VF_LOG_CRITICAL("Bad Alloc: {}", e.what()); - } catch (const std::exception &e) { - VF_LOG_CRITICAL("exception: {}", e.what()); - } catch (...) { - VF_LOG_CRITICAL("Unknown exception!"); - } + + try { + vf::logging::Logger::initializeLogger(); + vf::basics::ConfigurationFile config = + vf::basics::loadConfig(argc, argv, "./apps/gpu/DrivenCavityMultiGPU/drivencavity_1gpu.cfg"); + runVirtualFluids(config); + } catch (const std::exception& e) { + VF_LOG_WARNING("{}", e.what()); + return 1; } - MPI_Finalize(); return 0; } diff --git a/apps/gpu/SphereInChannel/SphereInChannel.cpp b/apps/gpu/SphereInChannel/SphereInChannel.cpp index bb63493de46e4df3e94e00036961a22072f14925..4bf4457c78dc6c608b5a224722e2e4032f8825ce 100644 --- a/apps/gpu/SphereInChannel/SphereInChannel.cpp +++ b/apps/gpu/SphereInChannel/SphereInChannel.cpp @@ -259,14 +259,9 @@ int main(int argc, char* argv[]) Simulation sim(para, cudaMemoryManager, communicator, *gridGenerator, &bcFactory, &scalingFactory); sim.run(); - } catch (const spdlog::spdlog_ex& ex) { - std::cout << "Log initialization failed: " << ex.what() << std::endl; - } catch (const std::bad_alloc& e) { - VF_LOG_CRITICAL("Bad Alloc: {}", e.what()); } catch (const std::exception& e) { - VF_LOG_CRITICAL("exception: {}", e.what()); - } catch (...) { - VF_LOG_CRITICAL("Unknown exception!"); + VF_LOG_WARNING("{}", e.what()); + return 1; } return 0; diff --git a/apps/gpu/SphereMultiGPU/SphereMultiGPU.cpp b/apps/gpu/SphereMultiGPU/SphereMultiGPU.cpp index f5a96c0373c87410f47c8c0dfeaaf3c03980772b..18df79ad745165a07755e65013bfa47b49f4b884 100755 --- a/apps/gpu/SphereMultiGPU/SphereMultiGPU.cpp +++ b/apps/gpu/SphereMultiGPU/SphereMultiGPU.cpp @@ -79,7 +79,7 @@ #include "gpu/core/Parameter/Parameter.h" #include "gpu/core/PreProcessor/PreProcessorFactory/PreProcessorFactoryImp.h" -void runVirtualFluids(const vf::basics::ConfigurationFile& config) +void run(const vf::basics::ConfigurationFile& config) { vf::parallel::Communicator& communicator = *vf::parallel::MPICommunicator::getInstance(); const auto numberOfProcesses = communicator.getNumberOfProcesses(); @@ -200,36 +200,21 @@ void runVirtualFluids(const vf::basics::ConfigurationFile& config) auto cudaMemoryManager = std::make_shared<CudaMemoryManager>(para); SPtr<GridProvider> gridGenerator = GridProvider::makeGridGenerator(gridBuilderFacade->getGridBuilder(), para, cudaMemoryManager, communicator); Simulation sim(para, cudaMemoryManager, communicator, *gridGenerator, &bcFactory, &scalingFactory); - - // run simulation sim.run(); } -int main(int argc, char *argv[]) +int main(int argc, char* argv[]) { - MPI_Init(&argc, &argv); - std::string str, str2, configFile; - - if (argv != NULL) { - - try { - VF_LOG_TRACE("For the default config path to work, execute the app from the project root."); - vf::basics::ConfigurationFile config = vf::basics::loadConfig(argc, argv, "./apps/gpu/SphereMultiGPU/sphere_1gpu.cfg"); - runVirtualFluids(config); - - ////////////////////////////////////////////////////////////////////////// - } catch (const spdlog::spdlog_ex &ex) { - std::cout << "Log initialization failed: " << ex.what() << std::endl; - } catch (const std::bad_alloc &e) { - VF_LOG_CRITICAL("Bad Alloc: {}", e.what()); - } catch (const std::exception &e) { - VF_LOG_CRITICAL("exception: {}", e.what()); - } catch (...) { - VF_LOG_CRITICAL("Unknown exception!"); - } + try { + vf::logging::Logger::initializeLogger(); + vf::basics::ConfigurationFile config = + vf::basics::loadConfig(argc, argv, "./apps/gpu/SphereMultiGPU/sphere_1gpu.cfg"); + run(config); + } catch (const std::exception& e) { + VF_LOG_WARNING("{}", e.what()); + return 1; } - MPI_Finalize(); return 0; } diff --git a/apps/gpu/TGV_3D/TGV_3D.cpp b/apps/gpu/TGV_3D/TGV_3D.cpp index 622add3afbc9225a18a029265f21ea3004d02dbb..84571931b720687b036d886b12de16660296cec9 100644 --- a/apps/gpu/TGV_3D/TGV_3D.cpp +++ b/apps/gpu/TGV_3D/TGV_3D.cpp @@ -272,72 +272,53 @@ void multipleLevel(const std::string& configPath) int main( int argc, char* argv[]) { - MPI_Init(&argc, &argv); std::string str, str2; - if ( argv != NULL ) - { - //str = static_cast<std::string>(argv[0]); - try - { - ////////////////////////////////////////////////////////////////////////// - std::string targetPath( __FILE__ ); + try { + ////////////////////////////////////////////////////////////////////////// + std::string targetPath(__FILE__); #ifdef _WIN32 - targetPath = targetPath.substr(0, targetPath.find_last_of('\\') + 1); + targetPath = targetPath.substr(0, targetPath.find_last_of('\\') + 1); #else - targetPath = targetPath.substr(0, targetPath.find_last_of('/') + 1); + targetPath = targetPath.substr(0, targetPath.find_last_of('/') + 1); #endif - ////////////////////////////////////////////////////////////////////////// - - if( cmdOptionExists( argv, argv+argc, "--Re" ) ) - reynoldsNumber = atof( getCmdOption( argv, argv+argc, "--Re" ) ); + ////////////////////////////////////////////////////////////////////////// - if( cmdOptionExists( argv, argv+argc, "--nx" ) ) - numberOfNodesX = atoi( getCmdOption( argv, argv+argc, "--nx" ) ); + if (cmdOptionExists(argv, argv + argc, "--Re")) + reynoldsNumber = atof(getCmdOption(argv, argv + argc, "--Re")); - if( cmdOptionExists( argv, argv+argc, "--dtPerL" ) ) - dtPerL = atoi( getCmdOption( argv, argv+argc, "--dtPerL" ) ); + if (cmdOptionExists(argv, argv + argc, "--nx")) + numberOfNodesX = atoi(getCmdOption(argv, argv + argc, "--nx")); - if( cmdOptionExists( argv, argv+argc, "--kernel" ) ) - kernel = getCmdOption( argv, argv+argc, "--kernel" ); + if (cmdOptionExists(argv, argv + argc, "--dtPerL")) + dtPerL = atoi(getCmdOption(argv, argv + argc, "--dtPerL")); - if( cmdOptionExists( argv, argv+argc, "--gpu" ) ) - gpuIndex = atoi( getCmdOption( argv, argv+argc, "--gpu" ) ); + if (cmdOptionExists(argv, argv + argc, "--kernel")) + kernel = getCmdOption(argv, argv + argc, "--kernel"); - if( cmdOptionExists( argv, argv+argc, "--useLimiter" ) ) - useLimiter = true; + if (cmdOptionExists(argv, argv + argc, "--gpu")) + gpuIndex = atoi(getCmdOption(argv, argv + argc, "--gpu")); - multipleLevel(targetPath + "tgv3d.cfg"); + if (cmdOptionExists(argv, argv + argc, "--useLimiter")) + useLimiter = true; - ////////////////////////////////////////////////////////////////////////// - } - catch (const std::bad_alloc& e) - { - std::cout << e.what() << std::flush; - //MPI_Abort(MPI_COMM_WORLD, -1); - } - catch (const std::exception& e) - { - std::cout << e.what() << std::flush; - //MPI_Abort(MPI_COMM_WORLD, -1); - } - catch (...) - { - std::cout << "unknown exeption" << std::endl; - } + multipleLevel(targetPath + "tgv3d.cfg"); - //std::cout << "\nConfiguration file must be set!: lbmgm <config file>" << std::endl << std::flush; - //MPI_Abort(MPI_COMM_WORLD, -1); + } catch (const std::exception& e) { + VF_LOG_WARNING("{}", e.what()); + return 1; } + // std::cout << "\nConfiguration file must be set!: lbmgm <config file>" << std::endl << std::flush; + // MPI_Abort(MPI_COMM_WORLD, -1); - /* - MPE_Init_log() & MPE_Finish_log() are NOT needed when - liblmpe.a is linked with this program. In that case, - MPI_Init() would have called MPE_Init_log() already. - */ + /* + MPE_Init_log() & MPE_Finish_log() are NOT needed when + liblmpe.a is linked with this program. In that case, + MPI_Init() would have called MPE_Init_log() already. + */ #if defined( MPI_LOGGING ) MPE_Init_log(); #endif @@ -351,7 +332,6 @@ int main( int argc, char* argv[]) MPE_Finish_log( "TestLog" ); #endif - MPI_Finalize(); return 0; } diff --git a/docs/pages/Build-and-Run.md b/docs/pages/Build-and-Run.md index e4c623acda14a0fe6257d2bdc7f2523185e70759..2196db39dd8a7e3083cb561f66f97d926494065c 100644 --- a/docs/pages/Build-and-Run.md +++ b/docs/pages/Build-and-Run.md @@ -3,7 +3,7 @@ # Build and Run -This guide describes how to start using and developing VirtualFluids in the terminal. Alternativly you can use Visual Studio to build an run it. +This guide describes how to start using and developing VirtualFluids in the terminal. Alternatively you can use Visual Studio to build an run it. ## Build @@ -58,9 +58,11 @@ For instance a simulation on the GPU containing a flow around a sphere can be st ``` +### Result files + The result files of this simulation are usually stored in: `./output/` -The result files of VirtualFluids are mostly in [VTK](https://kitware.github.io/vtk-examples/site/VTKFileFormats/) format. These files can be visualised with the free software [Paraview](https://www.paraview.org/). +The result files of VirtualFluids are mostly in [VTK](https://kitware.github.io/vtk-examples/site/VTKFileFormats/) format. These files can be visualized with the free software [Paraview](https://www.paraview.org/). -The CPU part generates a set of multiple output directories in the prescribed output path. The flow fields can be found in the _mq_ directory. To view the flow fields, it is most conveniant to open the _mq_collection.pvd_ file in Paraview. The _bc_ directory contains the boundary condition information, the _geo_ directory contains information on the geometry of the flow domain and the _blocks_ directory contains the block grid. +The CPU part generates a set of multiple output directories in the prescribed output path. The flow fields can be found in the _mq_ directory. To view the flow fields, it is most convenient to open the _mq_collection.pvd_ file in Paraview. The _bc_ directory contains the boundary condition information, the _geo_ directory contains information on the geometry of the flow domain and the _blocks_ directory contains the block grid. A GPU computation generates a the time series of output files directly in the output path. In Paraview these time series can be read directly. \ No newline at end of file diff --git a/src/basics/config/ConfigurationFile.cpp b/src/basics/config/ConfigurationFile.cpp index afd3ef3080467efc3b4ea4fc715fd362f76cbb75..81c657e7217edaa261461b08aaf2200d7dfc21a8 100644 --- a/src/basics/config/ConfigurationFile.cpp +++ b/src/basics/config/ConfigurationFile.cpp @@ -37,6 +37,7 @@ #include <fstream> #include <map> #include <sstream> +#include <stdexcept> #include <string> #include <vector> @@ -56,13 +57,17 @@ void ConfigurationFile::clear() data.clear(); } ////////////////////////////////////////////////////////////////////////// -bool ConfigurationFile::load(const std::string& file) +void ConfigurationFile::load(const std::string& file) { std::ifstream inFile(file.c_str()); if (!inFile.good()) { - UB_THROW(UbException(UB_EXARGS, "Cannot read configuration file " + file + "! Your current directory is " + - std::filesystem::current_path().string() + ".")); + const std::string error = "Cannot read configuration file " + file + "! Your current directory is " + + std::filesystem::current_path().string() + "\n" + + "For further information on how to run VirtualFluids please visit: " + "https://irmb.gitlab-pages.rz.tu-bs.de/VirtualFluids/build-and-run.html#run-the-examples"; + + throw std::invalid_argument(error); } while (inFile.good() && !inFile.eof()) { @@ -92,8 +97,6 @@ bool ConfigurationFile::load(const std::string& file) } } } - - return true; } ////////////////////////////////////////////////////////////////////////// diff --git a/src/basics/config/ConfigurationFile.h b/src/basics/config/ConfigurationFile.h index bf7293005053fe9cfb86f4acd8abfefe8a87cd2a..567a3fbe1710fd705836c6a6ff46441ca25b2029 100644 --- a/src/basics/config/ConfigurationFile.h +++ b/src/basics/config/ConfigurationFile.h @@ -98,7 +98,7 @@ public: void clear(); //! load a configuration file - bool load(const std::string& File); + void load(const std::string& File); //! check if value associated with given key exists bool contains(const std::string& key) const;