diff --git a/source/Applications/DLR-F16/F16BombadilTest10e-6.cfg b/source/Applications/DLR-F16/F16BombadilTest10e-6.cfg
index d94fc953bde419f8ad552c9f6d7aa608a66624ec..d9c64f2230860cdf5f68ea0b447ea2afe30a976c 100644
--- a/source/Applications/DLR-F16/F16BombadilTest10e-6.cfg
+++ b/source/Applications/DLR-F16/F16BombadilTest10e-6.cfg
@@ -9,7 +9,7 @@ zigZagTape = 2zackenbaender0.stl
numOfThreads = 4
availMem = 13e9
-refineLevel = 8 #10
+refineLevel = 1 #10
#blockNx = 7 8 8
#blockNx = 7 6 7
blockNx = 21 6 13
@@ -38,25 +38,26 @@ boundingBox = -0.90 1.20 0.035 0.065 -0.65 0.65
#deltaXfine = 0.000009765625
#deltaXfine = 0.005 #level 0
-#deltaXfine = 0.0025 #level 1
+deltaXfine = 0.0025 #level 1
#deltaXfine = 0.00125 #level 2
#deltaXfine = 0.000625 #level 3
#deltaXfine = 0.0003125 #level 4
#deltaXfine = 0.00015625 #level 5
#deltaXfine = 0.000078125 #level 6
#deltaXfine = 0.0000390625 #level 7
-deltaXfine = 0.00001953125 #level 8
+#deltaXfine = 0.00001953125 #level 8
#deltaXfine = 6.5e-6
startDistance = -1.0
refineDistance = 0.3
-restartStep = 500000
-restartStepStart = 500000
+newStart = false
+restartStep = 30
+restartStepStart = 20
-outTime = 1
-endTime = 15000
+outTime = 10
+endTime = 30
logToFile = false
diff --git a/source/Applications/DLR-F16/f16.cpp b/source/Applications/DLR-F16/f16.cpp
index adf9597c5dbfe440fa02360e5522cbf04b5158e2..02cc28c0a036a379d5da8d29972ae483dea38d6d 100644
--- a/source/Applications/DLR-F16/f16.cpp
+++ b/source/Applications/DLR-F16/f16.cpp
@@ -42,6 +42,7 @@ void run(string configname)
double refineDistance = config.getDouble("refineDistance");
double startDistance = config.getDouble("startDistance");
vector<double> nupsStep = config.getVector<double>("nupsStep");
+ bool newStart = config.getBool("newStart");
CommunicatorPtr comm = MPICommunicator::getInstance();
int myid = comm->getProcessID();
@@ -119,18 +120,6 @@ void run(string configname)
//Grid
//////////////////////////////////////////////////////////////////////////
Grid3DPtr grid(new Grid3D(comm));
- grid->setDeltaX(deltaXcoarse);
- grid->setBlockNX(blockNx[0], blockNx[1], blockNx[2]);
-
- GbObject3DPtr gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3));
- //gridCube->setCenterCoordinates(geo->getX1Centroid(), geo->getX2Centroid(), geo->getX3Centroid());
- if (myid == 0) GbSystem3D::writeGeoObject(gridCube.get(), pathOut + "/geo/gridCube", WbWriterVtkXmlASCII::getInstance());
- GenBlocksGridVisitor genBlocks(gridCube);
- grid->accept(genBlocks);
-
- grid->setPeriodicX1(false);
- grid->setPeriodicX2(true);
- grid->setPeriodicX3(false);
//BC adapters
BCAdapterPtr noSlipBCAdapter(new NoSlipBCAdapter());
@@ -163,13 +152,30 @@ void run(string configname)
//////////////////////////////////////////////////////////////////////////
//restart
- UbSchedulerPtr rSch(new UbScheduler(restartStep, restartStep));
- RestartCoProcessor rp(grid, rSch, comm, pathOut, RestartCoProcessor::TXT);
+ UbSchedulerPtr rSch(new UbScheduler(restartStep, restartStepStart));
+ //RestartCoProcessor rp(grid, rSch, comm, pathOut, RestartCoProcessor::TXT);
+ MPIIORestartCoProcessor rcp(grid, rSch, pathOut, comm);
//////////////////////////////////////////////////////////////////////////
- if (grid->getTimeStep() == 0)
+ //if (grid->getTimeStep() == 0)
+ if(newStart)
{
+ ////////////////////////////////////////////////////////////////////////
+ //define grid
+ //////////////////////////////////////////////////////////////////////////
+ grid->setDeltaX(deltaXcoarse);
+ grid->setBlockNX(blockNx[0], blockNx[1], blockNx[2]);
+
+ GbObject3DPtr gridCube(new GbCuboid3D(g_minX1, g_minX2, g_minX3, g_maxX1, g_maxX2, g_maxX3));
+ if (myid==0) GbSystem3D::writeGeoObject(gridCube.get(), pathOut+"/geo/gridCube", WbWriterVtkXmlASCII::getInstance());
+ GenBlocksGridVisitor genBlocks(gridCube);
+ grid->accept(genBlocks);
+
+ grid->setPeriodicX1(false);
+ grid->setPeriodicX2(true);
+ grid->setPeriodicX3(false);
+
if (myid == 0)
{
UBLOG(logINFO, "Parameters:");
@@ -520,7 +526,7 @@ void run(string configname)
//////////////////////////////////////
{
- WriteBlocksCoProcessor ppblocks(grid, UbSchedulerPtr(new UbScheduler(1)), pathOut, WbWriterVtkXmlBinary::getInstance(), comm);
+ WriteBlocksCoProcessor ppblocks(grid, UbSchedulerPtr(new UbScheduler(1)), pathOut, WbWriterVtkXmlASCII::getInstance(), comm);
ppblocks.process(2);
}
@@ -620,6 +626,14 @@ void run(string configname)
}
else
{
+ rcp.restart();
+ grid->setTimeStep(restartStepStart);
+
+ {
+ WriteBlocksCoProcessor ppblocks(grid, UbSchedulerPtr(new UbScheduler(1)), pathOut, WbWriterVtkXmlASCII::getInstance(), comm);
+ ppblocks.process(3);
+ }
+
InterpolationProcessorPtr iProcessor(new CompressibleOffsetInterpolationProcessor());
SetConnectorsBlockVisitor setConnsVisitor(comm, true, D3Q27System::ENDDIR, nuLB, iProcessor);
grid->accept(setConnsVisitor);
diff --git a/source/ThirdParty/Library/numerics/geometry3d/CoordinateTransformation3D.h b/source/ThirdParty/Library/numerics/geometry3d/CoordinateTransformation3D.h
index 9b0b83dec47590765f0d99e2cb37565d78b5cafd..e30495e61b3b7cf522ff2222d4dd78a3a631b85a 100644
--- a/source/ThirdParty/Library/numerics/geometry3d/CoordinateTransformation3D.h
+++ b/source/ThirdParty/Library/numerics/geometry3d/CoordinateTransformation3D.h
@@ -154,6 +154,8 @@ private:
bool active;
bool transformation;
+ friend class MPIIORestartCoProcessor;
+
friend class boost::serialization::access;
template<class Archive>
void serialize(Archive & ar, const unsigned int version)
diff --git a/source/VirtualFluids.h b/source/VirtualFluids.h
index 1a33e44b0e5ee971140cc0e40a91b6707d900228..ea027ec40d86432f3dec7ebbff7ffc6b1aafa88b 100644
--- a/source/VirtualFluids.h
+++ b/source/VirtualFluids.h
@@ -160,6 +160,7 @@
//#include <CoProcessors/MeanValuesCoProcessor.h>
#include <CoProcessors/TimeAveragedValuesCoProcessor.h>
#include <CoProcessors/InSituCatalystCoProcessor.h>
+#include <CoProcessors/MPIIORestartCoProcessor.h>
#include <LineTimeSeriesCoProcessor.h>
#include <IntegrateValuesHelper.h>
//#include <LBM/D3Q27CompactInterpolationProcessor.h>
diff --git a/source/VirtualFluidsCore/BoundaryConditions/BCArray3D.h b/source/VirtualFluidsCore/BoundaryConditions/BCArray3D.h
index 649584ae5f0f01dadc91d80fb136d97650ebf81b..2c44044244df36d43697d6305d469879e247f49e 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/BCArray3D.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/BCArray3D.h
@@ -99,6 +99,8 @@ private:
//////////////////////////////////////////////////////////////////////////
void deleteBC(std::size_t x1, std::size_t x2, std::size_t x3);
+ friend class MPIIORestartCoProcessor;
+
friend class boost::serialization::access;
template<class Archive>
void serialize(Archive & ar, const unsigned int version)
diff --git a/source/VirtualFluidsCore/BoundaryConditions/BoundaryConditions.h b/source/VirtualFluidsCore/BoundaryConditions/BoundaryConditions.h
index edbf155ceb314811ff90df2784afb286d967fa2a..3829c751cf10def0e0c09b00ebbe850ca0a91465 100644
--- a/source/VirtualFluidsCore/BoundaryConditions/BoundaryConditions.h
+++ b/source/VirtualFluidsCore/BoundaryConditions/BoundaryConditions.h
@@ -256,6 +256,8 @@ protected:
char algorithmType;
private:
+ friend class MPIIORestartCoProcessor;
+
friend class boost::serialization::access;
template<class Archive>
void serialize(Archive & ar, const unsigned int version)
diff --git a/source/VirtualFluidsCore/CoProcessors/MPIIORestartCoProcessor.cpp b/source/VirtualFluidsCore/CoProcessors/MPIIORestartCoProcessor.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..2d0bb5aefa3f2bf932c0350c29c5b09a9b64dae6
--- /dev/null
+++ b/source/VirtualFluidsCore/CoProcessors/MPIIORestartCoProcessor.cpp
@@ -0,0 +1,1036 @@
+#include "MPIIORestartCoProcessor.h"
+#include <boost/foreach.hpp>
+#include "D3Q27System.h"
+//#include "LBMKernel.h"
+#include "CompressibleCumulantLBMKernel.h"
+#include "D3Q27EsoTwist3DSplittedVector.h"
+
+//! \brief BLOCK_SIZE defines the quantity of the BoundaryCondition-structures written as one block to the file
+//! \details To avoid overflow in the parameter \a count of the function MPI_File_write_at
+//! structures BoundaryCondition are being written in blocks containing each of them BLOCK_SIZE structures
+#define BLOCK_SIZE 1024
+
+MPIIORestartCoProcessor::MPIIORestartCoProcessor(Grid3DPtr grid, UbSchedulerPtr s,
+ const std::string& path,
+ CommunicatorPtr comm) :
+ CoProcessor(grid, s),
+ path(path),
+ comm(comm)
+{
+ memset(&blockParamStr, 0, sizeof(blockParamStr));
+
+ //------------------------- define MPI types ---------------------------------
+
+ MPI_Datatype typesGP[3] = { MPI_DOUBLE, MPI_INT, MPI_CHAR };
+ int blocksGP[3] = { 34, 6, 5 };
+ MPI_Aint offsetsGP[3], lbGP, extentGP;
+
+ offsetsGP[0] = 0;
+ MPI_Type_get_extent(MPI_DOUBLE, &lbGP, &extentGP);
+ offsetsGP[1] = blocksGP[0] * extentGP;
+
+ MPI_Type_get_extent(MPI_INT, &lbGP, &extentGP);
+ offsetsGP[2] = offsetsGP[1] + blocksGP[1] * extentGP;
+
+ MPI_Type_struct(3, blocksGP, offsetsGP, typesGP, &gridParamType);
+ MPI_Type_commit(&gridParamType);
+
+ //-----------------------------------------------------------------------
+
+ MPI_Type_contiguous(41, MPI_INT, &blockParamType);
+ MPI_Type_commit(&blockParamType);
+
+ //-----------------------------------------------------------------------
+
+ MPI_Datatype types[3] = { MPI_DOUBLE, MPI_INT, MPI_CHAR };
+ int blocks[3] = { 2, 14, 3 };
+ MPI_Aint offsets[3], lb, extent;
+
+ offsets[0] = 0;
+ MPI_Type_get_extent(MPI_DOUBLE, &lb, &extent);
+ offsets[1] = blocks[0] * extent;
+
+ MPI_Type_get_extent(MPI_INT, &lb, &extent);
+ offsets[2] = offsets[1] + blocks[1] * extent;
+
+ MPI_Type_create_struct(3, blocks, offsets, types, &block3dType);
+ MPI_Type_commit(&block3dType);
+
+ //-----------------------------------------------------------------------
+
+ MPI_Type_contiguous(4, MPI_INT, &dataSetType);
+ MPI_Type_commit(&dataSetType);
+
+ //-----------------------------------------------------------------------
+
+ MPI_Datatype typesBC[3] = { MPI_LONG_LONG_INT, MPI_FLOAT, MPI_CHAR };
+ int blocksBC[3] = { 5, 38, 1 };
+ MPI_Aint offsetsBC[3], lbBC, extentBC;
+
+ offsetsBC[0] = 0;
+ MPI_Type_get_extent(MPI_LONG_LONG_INT, &lbBC, &extentBC);
+ offsetsBC[1] = blocksBC[0] * extentBC;
+
+ MPI_Type_get_extent(MPI_FLOAT, &lbBC, &extentBC);
+ offsetsBC[2] = blocksBC[1] * extentBC;
+
+ MPI_Type_struct(3, blocksBC, offsetsBC, typesBC, &boundCondType);
+ MPI_Type_commit(&boundCondType);
+
+ //---------------------------------------
+
+ MPI_Type_contiguous(BLOCK_SIZE, boundCondType, &boundCondType1000);
+ MPI_Type_commit(&boundCondType1000);
+
+ //---------------------------------------
+
+ MPI_Type_contiguous(6, MPI_INT, &boundCondTypeAdd);
+ MPI_Type_commit(&boundCondTypeAdd);
+
+}
+//////////////////////////////////////////////////////////////////////////
+MPIIORestartCoProcessor::~MPIIORestartCoProcessor()
+{
+ MPI_Type_free(&gridParamType);
+ MPI_Type_free(&blockParamType);
+ MPI_Type_free(&block3dType);
+ MPI_Type_free(&dataSetType);
+ MPI_Type_free(&dataSetDoubleType);
+ MPI_Type_free(&boundCondType);
+ MPI_Type_free(&boundCondType1000);
+ MPI_Type_free(&boundCondTypeAdd);
+ MPI_Type_free(&bcindexmatrixType);
+}
+
+//////////////////////////////////////////////////////////////////////////
+void MPIIORestartCoProcessor::process(double step)
+{
+ if(scheduler->isDue(step))
+ {
+ writeBlocks();
+ writeDataSet();
+ writeBoundaryConds();
+ }
+}
+//////////////////////////////////////////////////////////////////////////
+void MPIIORestartCoProcessor::writeBlocks()
+{
+ int rank, size;
+ MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+ MPI_Comm_size(MPI_COMM_WORLD, &size);
+ MPI_File file_handler;
+ int error = MPI_File_open(MPI_COMM_WORLD, "outputBlocks.bin", MPI_MODE_CREATE | MPI_MODE_RDWR , MPI_INFO_NULL, &file_handler);
+
+ int blocksCount = 0; // quantity of blocks in the grid, max 2147483648 blocks!
+ int minInitLevel = this->grid->getCoarsestInitializedLevel();
+ int maxInitLevel = this->grid->getFinestInitializedLevel();
+
+ std::vector<Block3DPtr> blocksVector[25]; // max 25 levels
+ for(int level = minInitLevel; level<=maxInitLevel; level++)
+ {
+ //grid->getBlocks(level, rank, blockVector[level]);
+ grid->getBlocks(level, blocksVector[level]);
+ blocksCount += static_cast<int>(blocksVector[level].size());
+ }
+
+ GridParam* gridParameters = new GridParam;
+ gridParameters->trafoParams[0] = grid->getCoordinateTransformator()->Tx1;
+ gridParameters->trafoParams[1] = grid->getCoordinateTransformator()->Tx2;
+ gridParameters->trafoParams[2] = grid->getCoordinateTransformator()->Tx3;
+ gridParameters->trafoParams[3] = grid->getCoordinateTransformator()->Sx1;
+ gridParameters->trafoParams[4] = grid->getCoordinateTransformator()->Sx2;
+ gridParameters->trafoParams[5] = grid->getCoordinateTransformator()->Sx3;
+ gridParameters->trafoParams[6] = grid->getCoordinateTransformator()->alpha;
+ gridParameters->trafoParams[7] = grid->getCoordinateTransformator()->beta;
+ gridParameters->trafoParams[8] = grid->getCoordinateTransformator()->gamma;
+
+ gridParameters->trafoParams[9] = grid->getCoordinateTransformator()->toX1factorX1;
+ gridParameters->trafoParams[10] = grid->getCoordinateTransformator()->toX1factorX2;
+ gridParameters->trafoParams[11] = grid->getCoordinateTransformator()->toX1factorX3;
+ gridParameters->trafoParams[12] = grid->getCoordinateTransformator()->toX1delta;
+ gridParameters->trafoParams[13] = grid->getCoordinateTransformator()->toX2factorX1;
+ gridParameters->trafoParams[14] = grid->getCoordinateTransformator()->toX2factorX2;
+ gridParameters->trafoParams[15] = grid->getCoordinateTransformator()->toX2factorX3;
+ gridParameters->trafoParams[16] = grid->getCoordinateTransformator()->toX2delta;
+ gridParameters->trafoParams[17] = grid->getCoordinateTransformator()->toX3factorX1;
+ gridParameters->trafoParams[18] = grid->getCoordinateTransformator()->toX3factorX2;
+ gridParameters->trafoParams[19] = grid->getCoordinateTransformator()->toX3factorX3;
+ gridParameters->trafoParams[20] = grid->getCoordinateTransformator()->toX3delta;
+
+ gridParameters->trafoParams[21] = grid->getCoordinateTransformator()->fromX1factorX1;
+ gridParameters->trafoParams[22] = grid->getCoordinateTransformator()->fromX1factorX2;
+ gridParameters->trafoParams[23] = grid->getCoordinateTransformator()->fromX1factorX3;
+ gridParameters->trafoParams[24] = grid->getCoordinateTransformator()->fromX1delta;
+ gridParameters->trafoParams[25] = grid->getCoordinateTransformator()->fromX2factorX1;
+ gridParameters->trafoParams[26] = grid->getCoordinateTransformator()->fromX2factorX2;
+ gridParameters->trafoParams[27] = grid->getCoordinateTransformator()->fromX2factorX3;
+ gridParameters->trafoParams[28] = grid->getCoordinateTransformator()->fromX2delta;
+ gridParameters->trafoParams[29] = grid->getCoordinateTransformator()->fromX3factorX1;
+ gridParameters->trafoParams[30] = grid->getCoordinateTransformator()->fromX3factorX2;
+ gridParameters->trafoParams[31] = grid->getCoordinateTransformator()->fromX3factorX3;
+ gridParameters->trafoParams[32] = grid->getCoordinateTransformator()->fromX3delta;
+
+ gridParameters->active = grid->getCoordinateTransformator()->active;
+ gridParameters->transformation = grid->getCoordinateTransformator()->transformation;
+
+ gridParameters->deltaX = grid->getDeltaX(minInitLevel);
+ UbTupleInt3 blocknx = grid->getBlockNX();
+ gridParameters->blockNx1 = val<1>(blocknx);
+ gridParameters->blockNx2 = val<2>(blocknx);
+ gridParameters->blockNx3 = val<3>(blocknx);
+ gridParameters->nx1 = grid->getNX1();
+ gridParameters->nx2 = grid->getNX2();
+ gridParameters->nx3 = grid->getNX3();
+ gridParameters->periodicX1 = grid->isPeriodicX1();
+ gridParameters->periodicX2 = grid->isPeriodicX2();
+ gridParameters->periodicX3 = grid->isPeriodicX3();
+
+//----------------------------------------------------------------------
+
+ Block3d* block3dArray = new Block3d[blocksCount];
+ bool firstBlock = true;
+ int ic = 0;
+ for(int level = minInitLevel; level <= maxInitLevel; level++)
+ {
+ BOOST_FOREACH(Block3DPtr block, blocksVector[level]) // all the blocks of the current level
+ {
+ if(firstBlock && block->getKernel()) // when first (any) valid block...
+ {
+ boost::shared_ptr< CbArray4D<LBMReal, IndexerX4X3X2X1> > AverageValuesArray3DPtr = block->getKernel()->getDataSet()->getAverageValues();
+ if (AverageValuesArray3DPtr)
+ {
+ blockParamStr.nx[0][0] = static_cast<int>(AverageValuesArray3DPtr->getNX1());
+ blockParamStr.nx[0][1] = static_cast<int>(AverageValuesArray3DPtr->getNX2());
+ blockParamStr.nx[0][2] = static_cast<int>(AverageValuesArray3DPtr->getNX3());
+ blockParamStr.nx[0][3] = static_cast<int>(AverageValuesArray3DPtr->getNX4());
+ }
+
+ boost::shared_ptr< CbArray4D<LBMReal, IndexerX4X3X2X1> > AverageVelocityArray3DPtr = block->getKernel()->getDataSet()->getAverageVelocity();
+ if (AverageVelocityArray3DPtr)
+ {
+ blockParamStr.nx[1][0] = static_cast<int>(AverageVelocityArray3DPtr->getNX1());
+ blockParamStr.nx[1][1] = static_cast<int>(AverageVelocityArray3DPtr->getNX2());
+ blockParamStr.nx[1][2] = static_cast<int>(AverageVelocityArray3DPtr->getNX3());
+ blockParamStr.nx[1][3] = static_cast<int>(AverageVelocityArray3DPtr->getNX4());
+ }
+
+ boost::shared_ptr< CbArray4D<LBMReal, IndexerX4X3X2X1> > AverageFluctArray3DPtr = block->getKernel()->getDataSet()->getAverageFluctuations();
+ if (AverageFluctArray3DPtr)
+ {
+ blockParamStr.nx[2][0] = static_cast<int>(AverageFluctArray3DPtr->getNX1());
+ blockParamStr.nx[2][1] = static_cast<int>(AverageFluctArray3DPtr->getNX2());
+ blockParamStr.nx[2][2] = static_cast<int>(AverageFluctArray3DPtr->getNX3());
+ blockParamStr.nx[2][3] = static_cast<int>(AverageFluctArray3DPtr->getNX4());
+ }
+
+ boost::shared_ptr< CbArray4D<LBMReal, IndexerX4X3X2X1> > AverageTripleArray3DPtr = block->getKernel()->getDataSet()->getAverageTriplecorrelations();
+ if (AverageTripleArray3DPtr)
+ {
+ blockParamStr.nx[3][0] = static_cast<int>(AverageTripleArray3DPtr->getNX1());
+ blockParamStr.nx[3][1] = static_cast<int>(AverageTripleArray3DPtr->getNX2());
+ blockParamStr.nx[3][2] = static_cast<int>(AverageTripleArray3DPtr->getNX3());
+ blockParamStr.nx[3][3] = static_cast<int>(AverageTripleArray3DPtr->getNX4());
+ }
+
+ boost::shared_ptr< CbArray4D<LBMReal, IndexerX4X3X2X1> > ShearStressValArray3DPtr = block->getKernel()->getDataSet()->getShearStressValues();
+ if (ShearStressValArray3DPtr)
+ {
+ blockParamStr.nx[4][0] = static_cast<int>(ShearStressValArray3DPtr->getNX1());
+ blockParamStr.nx[4][1] = static_cast<int>(ShearStressValArray3DPtr->getNX2());
+ blockParamStr.nx[4][2] = static_cast<int>(ShearStressValArray3DPtr->getNX3());
+ blockParamStr.nx[4][3] = static_cast<int>(ShearStressValArray3DPtr->getNX4());
+ }
+
+ boost::shared_ptr< CbArray3D<LBMReal, IndexerX3X2X1> > relaxationFactor3DPtr = block->getKernel()->getDataSet()->getRelaxationFactor();
+ if (relaxationFactor3DPtr)
+ {
+ blockParamStr.nx[5][0] = static_cast<int>(relaxationFactor3DPtr->getNX1());
+ blockParamStr.nx[5][1] = static_cast<int>(relaxationFactor3DPtr->getNX2());
+ blockParamStr.nx[5][2] = static_cast<int>(relaxationFactor3DPtr->getNX3());
+ blockParamStr.nx[5][3] = 1;
+ }
+
+ boost::shared_ptr< D3Q27EsoTwist3DSplittedVector > D3Q27EsoTwist3DSplittedVectorPtr = boost::dynamic_pointer_cast<D3Q27EsoTwist3DSplittedVector>(block->getKernel()->getDataSet()->getFdistributions());
+ CbArray4D<LBMReal, IndexerX4X3X2X1>::CbArray4DPtr localDistributions = D3Q27EsoTwist3DSplittedVectorPtr->getLocalDistributions();
+ if (localDistributions)
+ {
+ blockParamStr.nx[6][0] = static_cast<int>(localDistributions->getNX1());
+ blockParamStr.nx[6][1] = static_cast<int>(localDistributions->getNX2());
+ blockParamStr.nx[6][2] = static_cast<int>(localDistributions->getNX3());
+ blockParamStr.nx[6][3] = static_cast<int>(localDistributions->getNX4());
+ }
+
+ CbArray4D <LBMReal, IndexerX4X3X2X1>::CbArray4DPtr nonLocalDistributions = D3Q27EsoTwist3DSplittedVectorPtr->getNonLocalDistributions();
+ if (nonLocalDistributions)
+ {
+ blockParamStr.nx[7][0] = static_cast<int>(nonLocalDistributions->getNX1());
+ blockParamStr.nx[7][1] = static_cast<int>(nonLocalDistributions->getNX2());
+ blockParamStr.nx[7][2] = static_cast<int>(nonLocalDistributions->getNX3());
+ blockParamStr.nx[7][3] = static_cast<int>(nonLocalDistributions->getNX4());
+ }
+
+ CbArray3D<LBMReal, IndexerX3X2X1>::CbArray3DPtr zeroDistributions = D3Q27EsoTwist3DSplittedVectorPtr->getZeroDistributions();
+ if (zeroDistributions)
+ {
+ blockParamStr.nx[8][0] = static_cast<int>(zeroDistributions->getNX1());
+ blockParamStr.nx[8][1] = static_cast<int>(zeroDistributions->getNX2());
+ blockParamStr.nx[8][2] = static_cast<int>(zeroDistributions->getNX3());
+ blockParamStr.nx[8][3] = 1;
+ }
+
+ // ... than save some parameters that are equal in all blocks
+ blockParamStr.nx1 = static_cast<int>(block->getKernel()->getDataSet()->getFdistributions()->getNX1());
+ blockParamStr.nx2 = static_cast<int>(block->getKernel()->getDataSet()->getFdistributions()->getNX2());
+ blockParamStr.nx3 = static_cast<int>(block->getKernel()->getDataSet()->getFdistributions()->getNX3());
+
+ firstBlock = false;
+
+ // how many elements are in all arrays of DataSet (equal in all blocks)
+ int doubleCount = 0, temp;
+ for (int i = 0; i < 9; i++) // 9 arrays ( averageValues, averageVelocity, averageFluktuations,
+ { // averageTriplecorrelations, shearStressValues, relaxationFactor, 3 * fdistributions
+ temp = 1;
+ for (int ii = 0; ii < 4; ii++)
+ temp *= blockParamStr.nx[i][ii];
+ doubleCount += temp;
+ }
+ blockParamStr.doubleCountInBlock = doubleCount;
+
+ // the quantity of elements in the bcindexmatrix array (CbArray3D<int, IndexerX3X2X1>) in bcArray(BCArray3D) is always equal,
+ // this will be the size of the "write-read-block" in MPI_write_.../MPI_read-functions when writing/reading BoundConds
+ BCArray3D bcArr = block->getKernel()->getBCProcessor()->getBCArray();
+ blockParamStr.bcindexmatrix_count = static_cast<int>(bcArr.bcindexmatrix.getDataVector().size());
+ }
+
+ // save data describing the block
+ block3dArray[ic].x1 = block->getX1();
+ block3dArray[ic].x2 = block->getX2();
+ block3dArray[ic].x3 = block->getX3();
+ block3dArray[ic].bundle = block->getBundle();
+ block3dArray[ic].rank = block->getRank();
+ block3dArray[ic].lrank = block->getLocalRank();
+ block3dArray[ic].part = block->getPart();
+ block3dArray[ic].globalID = block->getGlobalID();
+ block3dArray[ic].localID = block->getLocalID();
+ block3dArray[ic].level = block->getLevel();
+ block3dArray[ic].interpolationFlagCF = block->getInterpolationFlagCF();
+ block3dArray[ic].interpolationFlagFC = block->getInterpolationFlagFC();
+ block3dArray[ic].counter = block->getMaxGlobalID();
+ block3dArray[ic].active = block->isActive();
+ if(block->getKernel())
+ {
+ block3dArray[ic].ghostLayerWidth = block->getKernel()->getGhostLayerWidth();
+ block3dArray[ic].collFactor = block->getKernel()->getCollisionFactor();
+ block3dArray[ic].deltaT = block->getKernel()->getDeltaT();
+ block3dArray[ic].compressible = block->getKernel()->getCompressible();
+ block3dArray[ic].withForcing = block->getKernel()->getWithForcing();
+ }
+ else
+ {
+ block3dArray[ic].ghostLayerWidth = 0;
+ block3dArray[ic].collFactor = 0.0;
+ block3dArray[ic].deltaT = 0.0;
+ block3dArray[ic].compressible = false;
+ block3dArray[ic].withForcing = false;
+ }
+
+ ic++;
+ }
+ }
+
+ // write to the file
+ // all processes calculate their offsets (quantity of bytes that the process is going to write)
+ // and notify the next process (with the rank = rank + 1)
+ size_t view_offset = size * sizeof(int);
+ size_t next_view_offset = 0;
+
+ if(size > 1)
+ {
+ if(rank == 0)
+ {
+ next_view_offset = view_offset + sizeof(GridParam) + sizeof(blockParam) + blocksCount * sizeof(Block3d);
+ MPI_Send(&next_view_offset, 1, MPI_LONG_LONG_INT, 1, 5, MPI_COMM_WORLD);
+ }
+ else
+ {
+ MPI_Recv(&view_offset, 1, MPI_LONG_LONG_INT, rank - 1, 5, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+ next_view_offset = view_offset + sizeof(GridParam) + sizeof(blockParam) + blocksCount * sizeof(Block3d);
+ if(rank < size - 1)
+ MPI_Send(&next_view_offset, 1, MPI_LONG_LONG_INT, rank + 1, 5, MPI_COMM_WORLD);
+ }
+ }
+
+ // each process writes the quantity of it's blocks
+ MPI_File_write_at(file_handler, rank * sizeof(int), &blocksCount, 1, MPI_INT, MPI_STATUS_IGNORE);
+ // each process writes parameters of the grid
+ MPI_File_write_at(file_handler, view_offset, gridParameters, 1, gridParamType, MPI_STATUS_IGNORE);
+ // each process writes common parameters of a block
+ MPI_File_write_at(file_handler, view_offset + sizeof(GridParam), &blockParamStr, 1, blockParamType, MPI_STATUS_IGNORE);
+ // each process writes it's blocks
+ MPI_File_write_at(file_handler, view_offset + sizeof(GridParam) + sizeof(blockParam), &block3dArray[0], blocksCount, block3dType, MPI_STATUS_IGNORE);
+ MPI_File_sync(file_handler);
+ MPI_File_close(&file_handler);
+
+ // register new MPI-types depending on the block-specific information
+ MPI_Type_contiguous(blockParamStr.doubleCountInBlock, MPI_DOUBLE, &dataSetDoubleType);
+ MPI_Type_commit(&dataSetDoubleType);
+
+ MPI_Type_contiguous(blockParamStr.bcindexmatrix_count, MPI_INT, &bcindexmatrixType);
+ MPI_Type_commit(&bcindexmatrixType);
+
+ delete[] block3dArray;
+ delete gridParameters;
+}
+
+void MPIIORestartCoProcessor::writeDataSet()
+{
+ int rank, size;
+ MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+ MPI_Comm_size(MPI_COMM_WORLD, &size);
+
+ int blocksCount = 0; // quantity of blocks in the grid, max 2147483648 blocks!
+
+ std::vector<Block3DPtr> blocksVector[25];
+ int minInitLevel = this->grid->getCoarsestInitializedLevel();
+ int maxInitLevel = this->grid->getFinestInitializedLevel();
+ for(int level = minInitLevel; level <= maxInitLevel;level++)
+ {
+ grid->getBlocks(level, rank, blocksVector[level]);
+ blocksCount += static_cast<int>(blocksVector[level].size());
+ }
+
+ dataSet* dataSetArray = new dataSet[blocksCount];
+ std::vector<double> doubleValuesArray; // double-values (arrays of f's) in all blocks
+
+ int ic = 0;
+ for(int level = minInitLevel; level<=maxInitLevel; level++)
+ {
+ BOOST_FOREACH(Block3DPtr block, blocksVector[level]) // blocks of the current level
+ {
+ dataSetArray[ic].x1 = block->getX1(); // coordinates of the block needed to find it while regenerating the grid
+ dataSetArray[ic].x2 = block->getX2();
+ dataSetArray[ic].x3 = block->getX3();
+ dataSetArray[ic].level = block->getLevel();
+
+ boost::shared_ptr< CbArray4D<LBMReal,IndexerX4X3X2X1> > AverageValuesArray3DPtr = block->getKernel()->getDataSet()->getAverageValues();
+ if(AverageValuesArray3DPtr && (blockParamStr.nx[0][0] > 0) && (blockParamStr.nx[0][1] > 0) && (blockParamStr.nx[0][2] > 0) && (blockParamStr.nx[0][3] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), AverageValuesArray3DPtr->getDataVector().begin(), AverageValuesArray3DPtr->getDataVector().end());
+
+ boost::shared_ptr< CbArray4D<LBMReal,IndexerX4X3X2X1> > AverageVelocityArray3DPtr = block->getKernel()->getDataSet()->getAverageVelocity();
+ if(AverageVelocityArray3DPtr && (blockParamStr.nx[1][0] > 0) && (blockParamStr.nx[1][1] > 0) && (blockParamStr.nx[1][2] > 0) && (blockParamStr.nx[1][3] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), AverageVelocityArray3DPtr->getDataVector().begin(), AverageVelocityArray3DPtr->getDataVector().end());
+
+ boost::shared_ptr< CbArray4D<LBMReal,IndexerX4X3X2X1> > AverageFluctArray3DPtr = block->getKernel()->getDataSet()->getAverageFluctuations();
+ if(AverageFluctArray3DPtr && (blockParamStr.nx[2][0] > 0) && (blockParamStr.nx[2][1] > 0) && (blockParamStr.nx[2][2] > 0) && (blockParamStr.nx[2][3] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), AverageFluctArray3DPtr->getDataVector().begin(), AverageFluctArray3DPtr->getDataVector().end());
+
+ boost::shared_ptr< CbArray4D<LBMReal,IndexerX4X3X2X1> > AverageTripleArray3DPtr = block->getKernel()->getDataSet()->getAverageTriplecorrelations();
+ if(AverageTripleArray3DPtr && (blockParamStr.nx[3][0] > 0) && (blockParamStr.nx[3][1] > 0) && (blockParamStr.nx[3][2] > 0) && (blockParamStr.nx[3][3] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), AverageTripleArray3DPtr->getDataVector().begin(), AverageTripleArray3DPtr->getDataVector().end());
+
+ boost::shared_ptr< CbArray4D<LBMReal,IndexerX4X3X2X1> > ShearStressValArray3DPtr = block->getKernel()->getDataSet()->getShearStressValues();
+ if(ShearStressValArray3DPtr && (blockParamStr.nx[4][0] > 0) && (blockParamStr.nx[4][1] > 0) && (blockParamStr.nx[4][2] > 0) && (blockParamStr.nx[4][3] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), ShearStressValArray3DPtr->getDataVector().begin(), ShearStressValArray3DPtr->getDataVector().end());
+
+ boost::shared_ptr< CbArray3D<LBMReal, IndexerX3X2X1> > RelaxationFactor3DPtr = block->getKernel()->getDataSet()->getRelaxationFactor();
+ if (RelaxationFactor3DPtr && (blockParamStr.nx[5][0] > 0) && (blockParamStr.nx[5][1] > 0) && (blockParamStr.nx[5][2] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), RelaxationFactor3DPtr->getDataVector().begin(), RelaxationFactor3DPtr->getDataVector().end());
+
+ boost::shared_ptr< D3Q27EsoTwist3DSplittedVector > D3Q27EsoTwist3DSplittedVectorPtr = boost::dynamic_pointer_cast<D3Q27EsoTwist3DSplittedVector>(block->getKernel()->getDataSet()->getFdistributions());
+ CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr localDistributions = D3Q27EsoTwist3DSplittedVectorPtr->getLocalDistributions();
+ if(localDistributions && (blockParamStr.nx[6][0] > 0) && (blockParamStr.nx[6][1] > 0) && (blockParamStr.nx[6][2] > 0) && (blockParamStr.nx[6][3] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), localDistributions->getDataVector().begin(), localDistributions->getDataVector().end());
+
+ CbArray4D <LBMReal,IndexerX4X3X2X1>::CbArray4DPtr nonLocalDistributions = D3Q27EsoTwist3DSplittedVectorPtr->getNonLocalDistributions();
+ if(nonLocalDistributions && (blockParamStr.nx[7][0] > 0) && (blockParamStr.nx[7][1] > 0) && (blockParamStr.nx[7][2] > 0) && (blockParamStr.nx[7][3] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), nonLocalDistributions->getDataVector().begin(), nonLocalDistributions->getDataVector().end());
+
+ CbArray3D<LBMReal,IndexerX3X2X1>::CbArray3DPtr zeroDistributions = D3Q27EsoTwist3DSplittedVectorPtr->getZeroDistributions();
+ if(zeroDistributions && (blockParamStr.nx[8][0] > 0) && (blockParamStr.nx[8][1] > 0) && (blockParamStr.nx[8][2] > 0))
+ doubleValuesArray.insert(doubleValuesArray.end(), zeroDistributions->getDataVector().begin(), zeroDistributions->getDataVector().end());
+
+ ic++;
+ }
+ }
+
+ // write to the file
+ // all processes calculate their offsets (quantity of bytes that the process is going to write)
+ // and notify the next process (with the rank = rank + 1)
+ size_t view_offset = size * sizeof(int);
+ size_t next_view_offset = 0;
+
+ if(size > 1)
+ {
+ if(rank == 0)
+ {
+ next_view_offset = view_offset + blocksCount * (sizeof(dataSet) + blockParamStr.doubleCountInBlock * sizeof(double));
+ MPI_Send(&next_view_offset, 1, MPI_LONG_LONG_INT, 1, 5, MPI_COMM_WORLD);
+ }
+ else
+ {
+ MPI_Recv(&view_offset, 1, MPI_LONG_LONG_INT, rank - 1, 5, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+ next_view_offset = view_offset + blocksCount * (sizeof(dataSet) + blockParamStr.doubleCountInBlock * sizeof(double));
+ if(rank < size - 1)
+ MPI_Send(&next_view_offset, 1, MPI_LONG_LONG_INT, rank + 1, 5, MPI_COMM_WORLD);
+ }
+ }
+
+ MPI_File file_handler;
+ int error = MPI_File_open(MPI_COMM_WORLD, "outputDataSet.bin", MPI_MODE_CREATE | MPI_MODE_RDWR , MPI_INFO_NULL, &file_handler);
+
+ // each process writes the quantity of it's blocks
+ MPI_File_write_at(file_handler, rank * sizeof(int), &blocksCount, 1, MPI_INT, MPI_STATUS_IGNORE);
+ // each process writes data identifying blocks
+ MPI_File_write_at(file_handler, view_offset, &dataSetArray[0], blocksCount, dataSetType, MPI_STATUS_IGNORE);
+ // each process writes the dataSet arrays
+ MPI_File_write_at(file_handler, view_offset + blocksCount * sizeof(dataSet), &doubleValuesArray[0], blocksCount, dataSetDoubleType, MPI_STATUS_IGNORE);
+ MPI_File_sync(file_handler);
+ MPI_File_close(&file_handler);
+
+ delete[] dataSetArray;
+}
+
+void MPIIORestartCoProcessor::writeBoundaryConds()
+{
+ int rank, size;
+ MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+ MPI_Comm_size(MPI_COMM_WORLD, &size);
+
+ int blocksCount = 0; // quantity of blocks in the grid, max 2147483648 blocks!
+ size_t count_boundCond = 0; // how many BoundaryConditions in all blocks
+ int count_indexContainer = 0; // how many indexContainer-values in all blocks
+ size_t byteCount = 0; // how many bytes writes this process in the file
+
+ std::vector<Block3DPtr> blocksVector[25];
+ int minInitLevel = this->grid->getCoarsestInitializedLevel();
+ int maxInitLevel = this->grid->getFinestInitializedLevel();
+ for(int level = minInitLevel; level <= maxInitLevel;level++)
+ {
+ grid->getBlocks(level, rank, blocksVector[level]);
+ blocksCount += static_cast<int>(blocksVector[level].size());
+ }
+
+ BCAdd* bcAddArray = new BCAdd[blocksCount];
+ std::vector<BoundaryCondition> bcVector;
+ std::vector<int> bcindexmatrixV;
+ std::vector<int> indexContainerV;
+
+ int ic = 0;
+ for(int level = minInitLevel; level <= maxInitLevel; level++)
+ {
+ BOOST_FOREACH(Block3DPtr block, blocksVector[level]) // all the blocks of the current level
+ {
+ BCArray3D bcArr = block->getKernel()->getBCProcessor()->getBCArray();
+
+ bcAddArray[ic].x1 = block->getX1(); // coordinates of the block needed to find it while regenerating the grid
+ bcAddArray[ic].x2 = block->getX2();
+ bcAddArray[ic].x3 = block->getX3();
+ bcAddArray[ic].level = block->getLevel();
+ bcAddArray[ic].boundCond_count = 0; // how many BoundaryConditions in this block
+ bcAddArray[ic].indexContainer_count = 0; // how many indexContainer-values in this block
+
+ for(int bc = 0; bc < bcArr.getBCVectorSize(); bc++)
+ {
+ BoundaryCondition* bouCond = new BoundaryCondition();
+ if(bcArr.bcvector[bc] == NULL)
+ {
+ memset(bouCond, 0, sizeof(BoundaryCondition));
+ }
+ else
+ {
+ bouCond->noslipBoundaryFlags = bcArr.bcvector[bc]->getNoSlipBoundary();
+ bouCond->slipBoundaryFlags = bcArr.bcvector[bc]->getSlipBoundary();
+ bouCond->velocityBoundaryFlags = bcArr.bcvector[bc]->getVelocityBoundary();
+ bouCond->densityBoundaryFlags = bcArr.bcvector[bc]->getDensityBoundary();
+ bouCond->wallModelBoundaryFlags = bcArr.bcvector[bc]->getWallModelBoundary();
+ bouCond->bcVelocityX1 = bcArr.bcvector[bc]->getBoundaryVelocityX1();
+ bouCond->bcVelocityX2 = bcArr.bcvector[bc]->getBoundaryVelocityX2();
+ bouCond->bcVelocityX3 = bcArr.bcvector[bc]->getBoundaryVelocityX3();
+ bouCond->bcDensity = bcArr.bcvector[bc]->getBoundaryDensity();
+ bouCond->bcLodiDensity = bcArr.bcvector[bc]->getDensityLodiDensity();
+ bouCond->bcLodiVelocityX1 = bcArr.bcvector[bc]->getDensityLodiVelocityX1();
+ bouCond->bcLodiVelocityX2 = bcArr.bcvector[bc]->getDensityLodiVelocityX2();
+ bouCond->bcLodiVelocityX3 = bcArr.bcvector[bc]->getDensityLodiVelocityX3();
+ bouCond->bcLodiLentgh = bcArr.bcvector[bc]->getDensityLodiLength();
+ bouCond->nx1 = bcArr.bcvector[bc]->nx1;
+ bouCond->nx2 = bcArr.bcvector[bc]->nx2;
+ bouCond->nx3 = bcArr.bcvector[bc]->nx3;
+ for(int iq = 0; iq < 26; iq++)
+ bouCond->q[iq] = bcArr.bcvector[bc]->getQ(iq);
+ bouCond->algorithmType = bcArr.bcvector[bc]->getBcAlgorithmType();
+ }
+
+ bcVector.push_back(*bouCond);
+ bcAddArray[ic].boundCond_count++;
+ count_boundCond++;
+ }
+
+ bcindexmatrixV.insert(bcindexmatrixV.end(), bcArr.bcindexmatrix.getDataVector().begin(), bcArr.bcindexmatrix.getDataVector().end());
+
+ indexContainerV.insert(indexContainerV.end(), bcArr.bcindexmatrix.getDataVector().begin(), bcArr.bcindexmatrix.getDataVector().end());
+ bcAddArray[ic].indexContainer_count = static_cast<int>(bcArr.indexContainer.size());
+ count_indexContainer += bcAddArray[ic].indexContainer_count;
+
+ ic++;
+ }
+ }
+
+ //how many "big blocks" of BLOCK_SIZE size can by formed
+ int bcBlockCount = (int)(count_boundCond / BLOCK_SIZE);
+ if(bcBlockCount * BLOCK_SIZE < count_boundCond)
+ bcBlockCount += 1;
+ for (int i = (int)count_boundCond; i < bcBlockCount * BLOCK_SIZE; i++)
+ {
+ BoundaryCondition* bouCond = new BoundaryCondition();
+ memset(bouCond, 0, sizeof(BoundaryCondition));
+ bcVector.push_back(*bouCond);
+ }
+
+ byteCount = bcBlockCount * BLOCK_SIZE * sizeof(BoundaryCondition) + blocksCount * sizeof(BCAdd) + sizeof(int) * (blocksCount * blockParamStr.bcindexmatrix_count + count_indexContainer);
+
+ // write to the file
+ // all processes calculate their offsets (quantity of bytes that the process is going to write)
+ // and notify the next process (with the rank = rank + 1)
+ size_t view_offset = 3 * size * sizeof(int);
+ size_t next_view_offset = 0;
+
+ if(size > 1)
+ {
+ if(rank == 0)
+ {
+ next_view_offset = view_offset + byteCount;
+ MPI_Send(&next_view_offset, 1, MPI_LONG_LONG_INT, 1, 5, MPI_COMM_WORLD);
+ }
+ else
+ {
+ MPI_Recv(&view_offset, 1, MPI_LONG_LONG_INT, rank - 1, 5, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+ next_view_offset = view_offset + byteCount;
+ if(rank < size - 1)
+ MPI_Send(&next_view_offset, 1, MPI_LONG_LONG_INT, rank + 1, 5, MPI_COMM_WORLD);
+ }
+ }
+
+ MPI_File file_handler;
+ int error = MPI_File_open(MPI_COMM_WORLD, "outputBoundCond.bin", MPI_MODE_CREATE | MPI_MODE_RDWR , MPI_INFO_NULL, &file_handler);
+
+ // each process writes the quantity of it's blocks
+ MPI_File_write_at(file_handler, rank * sizeof(int), &blocksCount, 1, MPI_INT, MPI_STATUS_IGNORE); // blocks quantity
+ // each process writes the quantity of "big blocks" of BLOCK_SIZE of boundary conditions
+ MPI_File_write_at(file_handler, (rank + size) * sizeof(int), &bcBlockCount, 1, MPI_INT, MPI_STATUS_IGNORE); // quantity of BoundConds / BLOCK_SIZE
+ // each process writes the quantity of indexContainer elements in all blocks
+ MPI_File_write_at(file_handler, (rank + 2 * size) * sizeof(int), &count_indexContainer, 1, MPI_INT, MPI_STATUS_IGNORE); // quantity of indexContainer
+
+ // each process writes data identifying the blocks
+ MPI_File_write_at(file_handler, view_offset, &bcAddArray[0], blocksCount, boundCondTypeAdd, MPI_STATUS_IGNORE);
+ // each process writes boundary conditions
+ MPI_File_write_at(file_handler, view_offset + blocksCount * sizeof(BCAdd), &bcVector[0], bcBlockCount, boundCondType1000, MPI_STATUS_IGNORE);
+ // each process writes bcindexmatrix values
+ MPI_File_write_at(file_handler, view_offset + blocksCount * sizeof(BCAdd) + bcBlockCount * BLOCK_SIZE * sizeof(BoundaryCondition), &bcindexmatrixV[0], blocksCount, bcindexmatrixType, MPI_STATUS_IGNORE);
+ // each process writes indexContainer values
+ MPI_File_write_at(file_handler, view_offset + blocksCount * sizeof(BCAdd) + bcBlockCount * BLOCK_SIZE * sizeof(BoundaryCondition) + blocksCount * blockParamStr.bcindexmatrix_count * sizeof(int), &indexContainerV[0], count_indexContainer, MPI_INT, MPI_STATUS_IGNORE);
+ MPI_File_sync(file_handler);
+ MPI_File_close(&file_handler);
+
+ delete [] bcAddArray;
+}
+
+//------------------------------------------- READ -----------------------------------------------
+void MPIIORestartCoProcessor::restart()
+{
+ readBlocks();
+ readDataSet();
+ readBoundaryConds();
+
+ this->reconnect(grid);
+}
+
+void MPIIORestartCoProcessor::readBlocks()
+{
+ int rank, size;
+ MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+ MPI_Comm_size(MPI_COMM_WORLD, &size);
+
+ MPI_File file_handler;
+ int error = MPI_File_open(MPI_COMM_WORLD, "outputBlocks.bin", MPI_MODE_CREATE | MPI_MODE_RDWR , MPI_INFO_NULL, &file_handler);
+
+ // read count of blocks
+ int blocksCount = 0;
+ MPI_File_read_at(file_handler, rank * sizeof(int), &blocksCount, 1, MPI_INT, MPI_STATUS_IGNORE);
+ Block3d* block3dArray = new Block3d[blocksCount];
+
+ // calculate the read offset
+ size_t read_offset = size * sizeof(int);
+ size_t next_read_offset = 0;
+
+ if(size > 1)
+ {
+ if(rank == 0)
+ {
+ next_read_offset = read_offset + sizeof(GridParam) + sizeof(blockParam) + blocksCount * sizeof(Block3d);
+ MPI_Send(&next_read_offset, 1, MPI_LONG_LONG_INT, 1, 5, MPI_COMM_WORLD);
+ }
+ else
+ {
+ MPI_Recv(&read_offset, 1, MPI_LONG_LONG_INT, rank - 1, 5, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+ next_read_offset = read_offset + sizeof(GridParam) + sizeof(blockParam) + blocksCount * sizeof(Block3d);
+ if(rank < size - 1)
+ MPI_Send(&next_read_offset, 1, MPI_LONG_LONG_INT, rank + 1, 5, MPI_COMM_WORLD);
+ }
+ }
+
+ GridParam* gridParameters = new GridParam;
+
+ // read parameters of the grid
+ MPI_File_read_at(file_handler, read_offset, gridParameters, 1, gridParamType, MPI_STATUS_IGNORE);
+ // read parameters of a block
+ MPI_File_read_at(file_handler, read_offset + sizeof(GridParam), &blockParamStr, 1, blockParamType, MPI_STATUS_IGNORE);
+ // read all the blocks
+ MPI_File_read_at(file_handler, read_offset + sizeof(GridParam) + sizeof(blockParam), &block3dArray[0], blocksCount, block3dType, MPI_STATUS_IGNORE);
+ MPI_File_sync(file_handler);
+
+ MPI_File_close(&file_handler);
+
+ // clear the grid
+ std::vector<Block3DPtr> blocksVector;
+ grid->getBlocks(0, blocksVector);
+ int del = 0;
+ BOOST_FOREACH(Block3DPtr block, blocksVector)
+ {
+ grid->deleteBlock(block);
+ del++;
+ }
+
+ // restore the grid
+ CoordinateTransformation3DPtr trafo(new CoordinateTransformation3D());
+ trafo->Tx1 = gridParameters->trafoParams[0];
+ trafo->Tx2 = gridParameters->trafoParams[1];
+ trafo->Tx3 = gridParameters->trafoParams[2];
+ trafo->Sx1 = gridParameters->trafoParams[3];
+ trafo->Sx2 = gridParameters->trafoParams[4];
+ trafo->Sx3 = gridParameters->trafoParams[5];
+ trafo->alpha = gridParameters->trafoParams[6];
+ trafo->beta = gridParameters->trafoParams[7];
+ trafo->gamma = gridParameters->trafoParams[8];
+
+ trafo->toX1factorX1 = gridParameters->trafoParams[9];
+ trafo->toX1factorX2 = gridParameters->trafoParams[10];
+ trafo->toX1factorX3 = gridParameters->trafoParams[11];
+ trafo->toX1delta = gridParameters->trafoParams[12];
+ trafo->toX2factorX1 = gridParameters->trafoParams[13];
+ trafo->toX2factorX2 = gridParameters->trafoParams[14];
+ trafo->toX2factorX3 = gridParameters->trafoParams[15];
+ trafo->toX2delta = gridParameters->trafoParams[16];
+ trafo->toX3factorX1 = gridParameters->trafoParams[17];
+ trafo->toX3factorX2 = gridParameters->trafoParams[18];
+ trafo->toX3factorX3 = gridParameters->trafoParams[19];
+ trafo->toX3delta = gridParameters->trafoParams[20];
+
+ trafo->fromX1factorX1 = gridParameters->trafoParams[21];
+ trafo->fromX1factorX2 = gridParameters->trafoParams[22];
+ trafo->fromX1factorX3 = gridParameters->trafoParams[23];
+ trafo->fromX1delta = gridParameters->trafoParams[24];
+ trafo->fromX2factorX1 = gridParameters->trafoParams[25];
+ trafo->fromX2factorX2 = gridParameters->trafoParams[26];
+ trafo->fromX2factorX3 = gridParameters->trafoParams[27];
+ trafo->fromX2delta = gridParameters->trafoParams[28];
+ trafo->fromX3factorX1 = gridParameters->trafoParams[29];
+ trafo->fromX3factorX2 = gridParameters->trafoParams[30];
+ trafo->fromX3factorX3 = gridParameters->trafoParams[31];
+ trafo->fromX3delta = gridParameters->trafoParams[32];
+
+ trafo->active = gridParameters->active;
+ trafo->transformation = gridParameters->transformation;
+
+ grid->setCoordinateTransformator(trafo);
+
+ grid->setDeltaX(gridParameters->deltaX);
+ grid->setBlockNX(gridParameters->blockNx1, gridParameters->blockNx2, gridParameters->blockNx3);
+ grid->setNX1(gridParameters->nx1);
+ grid->setNX2(gridParameters->nx2);
+ grid->setNX3(gridParameters->nx3);
+ grid->setPeriodicX1(gridParameters->periodicX1);
+ grid->setPeriodicX2(gridParameters->periodicX2);
+ grid->setPeriodicX3(gridParameters->periodicX3);
+
+ // regenerate blocks
+ for(int n = 0; n < blocksCount; n++)
+ {
+ LBMKernelPtr kernel(new CompressibleCumulantLBMKernel(block3dArray[n].x1, block3dArray[n].x2, block3dArray[n].x3, CompressibleCumulantLBMKernel::NORMAL));
+ kernel->setCollisionFactor(block3dArray[n].collFactor);
+ kernel->setGhostLayerWidth(block3dArray[n].ghostLayerWidth);
+ kernel->setCompressible(block3dArray[n].compressible);
+ kernel->setWithForcing(block3dArray[n].withForcing);
+ kernel->setDeltaT(block3dArray[n].deltaT);
+
+ Block3DPtr block(new Block3D(block3dArray[n].x1, block3dArray[n].x2, block3dArray[n].x3, block3dArray[n].level));
+ block->setActive(block3dArray[n].active);
+ block->setBundle(block3dArray[n].bundle);
+ block->setRank(block3dArray[n].rank);
+ block->setLocalRank(block3dArray[n].lrank);
+ block->setGlobalID(block3dArray[n].globalID);
+ block->setLocalID(block3dArray[n].localID);
+ block->setPart(block3dArray[n].part),
+ block->setLevel(block3dArray[n].level);
+ //block->setInterpolationFlagCF(block3dArray[n].interpolationFlagCF);
+ //block->setInterpolationFlagFC(block3dArray[n].interpolationFlagFC);
+ block->interpolationFlagCF = block3dArray[n].interpolationFlagCF;
+ block->interpolationFlagFC = block3dArray[n].interpolationFlagFC;
+
+ if(rank == block3dArray[n].rank)
+ block->setKernel(kernel);
+
+ grid->addBlock(block);
+ }
+
+ // define MPI_types depending on the block-specific information
+ MPI_Type_contiguous(blockParamStr.doubleCountInBlock, MPI_DOUBLE, &dataSetDoubleType);
+ MPI_Type_commit(&dataSetDoubleType);
+
+ MPI_Type_contiguous(blockParamStr.bcindexmatrix_count, MPI_INT, &bcindexmatrixType);
+ MPI_Type_commit(&bcindexmatrixType);
+
+ delete gridParameters;
+ delete [] block3dArray;
+}
+
+void MPIIORestartCoProcessor::readDataSet()
+{
+ int rank, size;
+ MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+ MPI_Comm_size(MPI_COMM_WORLD, &size);
+
+ MPI_File file_handler;
+ int error = MPI_File_open(MPI_COMM_WORLD, "outputDataSet.bin", MPI_MODE_CREATE | MPI_MODE_RDWR , MPI_INFO_NULL, &file_handler);
+
+ // read count of blocks
+ int blocksCount = 0;
+ MPI_File_read_at(file_handler, rank * sizeof(int), &blocksCount, 1, MPI_INT, MPI_STATUS_IGNORE);
+
+ dataSet* dataSetArray = new dataSet[blocksCount];
+ std::vector<double> doubleValuesArray(blocksCount * blockParamStr.doubleCountInBlock); // double-values in all blocks
+
+ // calculate the read offset
+ size_t read_offset = size * sizeof(int);
+ size_t next_read_offset = 0;
+
+ if(size > 1)
+ {
+ if(rank == 0)
+ {
+ next_read_offset = read_offset + blocksCount * (sizeof(dataSet) + blockParamStr.doubleCountInBlock * sizeof(double));
+ MPI_Send(&next_read_offset, 1, MPI_LONG_LONG_INT, 1, 5, MPI_COMM_WORLD);
+ }
+ else
+ {
+ MPI_Recv(&read_offset, 1, MPI_LONG_LONG_INT, rank - 1, 5, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+ next_read_offset = read_offset + blocksCount * (sizeof(dataSet) + blockParamStr.doubleCountInBlock * sizeof(double));
+ if(rank < size - 1)
+ MPI_Send(&next_read_offset, 1, MPI_LONG_LONG_INT, rank + 1, 5, MPI_COMM_WORLD);
+ }
+ }
+
+ MPI_File_read_at(file_handler, read_offset, &dataSetArray[0], blocksCount, dataSetType, MPI_STATUS_IGNORE);
+ MPI_File_read_at(file_handler, read_offset + blocksCount * sizeof(dataSet), &doubleValuesArray[0], blocksCount, dataSetDoubleType, MPI_STATUS_IGNORE);
+ MPI_File_sync(file_handler);
+ MPI_File_close(&file_handler);
+
+ size_t index = 0, nextVectorSize = 0;
+ std::vector<double> vectorsOfValues[9];
+ for (int n = 0; n < blocksCount; n++)
+ {
+ for(int b = 0; b < 9; b++) // assign approciate vectors for 9 dataSet arrays
+ {
+ nextVectorSize = blockParamStr.nx[b][0] * blockParamStr.nx[b][1] * blockParamStr.nx[b][2] * blockParamStr.nx[b][3];
+ vectorsOfValues[b].assign(doubleValuesArray.data() + index, doubleValuesArray.data() + index + nextVectorSize);
+ index += nextVectorSize;
+ }
+
+ // fill dataSet arrays
+ AverageValuesArray3DPtr mAverageValues;
+ if((blockParamStr.nx[0][0] == 0)&&(blockParamStr.nx[0][1] == 0)&&(blockParamStr.nx[0][2] == 0)&&(blockParamStr.nx[0][3] == 0))
+ mAverageValues = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr();
+ else
+ mAverageValues = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr(new CbArray4D<LBMReal,IndexerX4X3X2X1>(vectorsOfValues[0], blockParamStr.nx[0][0], blockParamStr.nx[0][1], blockParamStr.nx[0][2], blockParamStr.nx[0][3]));
+
+ AverageValuesArray3DPtr mAverageVelocity;
+ if((blockParamStr.nx[1][0] == 0)&&(blockParamStr.nx[1][1] == 0)&&(blockParamStr.nx[1][2] == 0)&&(blockParamStr.nx[1][3] == 0))
+ mAverageVelocity = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr();
+ else
+ mAverageVelocity = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr(new CbArray4D<LBMReal,IndexerX4X3X2X1>(vectorsOfValues[1], blockParamStr.nx[1][0], blockParamStr.nx[1][1], blockParamStr.nx[1][2], blockParamStr.nx[1][3]));
+
+ AverageValuesArray3DPtr mAverageFluktuations;
+ if((blockParamStr.nx[2][0] == 0)&&(blockParamStr.nx[2][1] == 0)&&(blockParamStr.nx[2][2] == 0)&&(blockParamStr.nx[2][3] == 0))
+ mAverageFluktuations = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr();
+ else
+ mAverageFluktuations = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr(new CbArray4D<LBMReal,IndexerX4X3X2X1>(vectorsOfValues[2], blockParamStr.nx[2][0], blockParamStr.nx[2][1], blockParamStr.nx[2][2], blockParamStr.nx[2][3]));
+
+ AverageValuesArray3DPtr mAverageTriplecorrelations;
+ if((blockParamStr.nx[3][0] == 0)&&(blockParamStr.nx[3][1] == 0)&&(blockParamStr.nx[3][2] == 0)&&(blockParamStr.nx[3][3] == 0))
+ mAverageTriplecorrelations = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr();
+ else
+ mAverageTriplecorrelations = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr(new CbArray4D<LBMReal,IndexerX4X3X2X1>(vectorsOfValues[3], blockParamStr.nx[3][0], blockParamStr.nx[3][1], blockParamStr.nx[3][2], blockParamStr.nx[3][3]));
+
+ ShearStressValuesArray3DPtr mShearStressValues;
+ if((blockParamStr.nx[4][0] == 0)&&(blockParamStr.nx[4][1] == 0)&&(blockParamStr.nx[4][2] == 0)&&(blockParamStr.nx[4][3] == 0))
+ mShearStressValues = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr();
+ else
+ mShearStressValues = CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr(new CbArray4D<LBMReal,IndexerX4X3X2X1>(vectorsOfValues[4], blockParamStr.nx[4][0], blockParamStr.nx[4][1], blockParamStr.nx[4][2], blockParamStr.nx[4][3]));
+
+ RelaxationFactorArray3DPtr mRelaxationFactor;
+ if ((blockParamStr.nx[5][0] == 0) && (blockParamStr.nx[5][1] == 0) && (blockParamStr.nx[5][2] == 0))
+ mRelaxationFactor = CbArray3D<LBMReal, IndexerX3X2X1>::CbArray3DPtr();
+ else
+ mRelaxationFactor = CbArray3D<LBMReal, IndexerX3X2X1>::CbArray3DPtr(new CbArray3D<LBMReal, IndexerX3X2X1>(vectorsOfValues[5], blockParamStr.nx[5][0], blockParamStr.nx[5][1], blockParamStr.nx[5][2]));
+
+ DistributionArray3DPtr mFdistributions(new D3Q27EsoTwist3DSplittedVector(blockParamStr.nx1, blockParamStr.nx2, blockParamStr.nx3, -999.0));
+
+ boost::dynamic_pointer_cast<D3Q27EsoTwist3DSplittedVector>(mFdistributions)->setLocalDistributions( CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr(new CbArray4D<LBMReal,IndexerX4X3X2X1>(vectorsOfValues[6], blockParamStr.nx[6][0], blockParamStr.nx[6][1], blockParamStr.nx[6][2], blockParamStr.nx[6][3])));
+ boost::dynamic_pointer_cast<D3Q27EsoTwist3DSplittedVector>(mFdistributions)->setNonLocalDistributions( CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr(new CbArray4D<LBMReal,IndexerX4X3X2X1>(vectorsOfValues[7], blockParamStr.nx[7][0], blockParamStr.nx[7][1], blockParamStr.nx[7][2], blockParamStr.nx[7][3])));
+ boost::dynamic_pointer_cast<D3Q27EsoTwist3DSplittedVector>(mFdistributions)->setZeroDistributions( CbArray3D<LBMReal,IndexerX3X2X1>::CbArray3DPtr(new CbArray3D<LBMReal,IndexerX3X2X1>(vectorsOfValues[8], blockParamStr.nx[8][0], blockParamStr.nx[8][1], blockParamStr.nx[8][2])));
+
+ boost::dynamic_pointer_cast<D3Q27EsoTwist3DSplittedVector>(mFdistributions)->setNX1(blockParamStr.nx1);
+ boost::dynamic_pointer_cast<D3Q27EsoTwist3DSplittedVector>(mFdistributions)->setNX2(blockParamStr.nx2);
+ boost::dynamic_pointer_cast<D3Q27EsoTwist3DSplittedVector>(mFdistributions)->setNX3(blockParamStr.nx3);
+
+ DataSet3DPtr dataSetPtr = DataSet3DPtr(new DataSet3D());
+ dataSetPtr->setAverageValues(mAverageValues);
+ dataSetPtr->setAverageVelocity(mAverageVelocity);
+ dataSetPtr->setAverageFluctuations(mAverageFluktuations);
+ dataSetPtr->setAverageTriplecorrelations(mAverageTriplecorrelations);
+ dataSetPtr->setShearStressValues(mShearStressValues);
+ dataSetPtr->setRelaxationFactor(mRelaxationFactor);
+ dataSetPtr->setFdistributions(mFdistributions);
+
+ // find the nesessary block and fill it
+ Block3DPtr block = grid->getBlock(dataSetArray[n].x1, dataSetArray[n].x2, dataSetArray[n].x3, dataSetArray[n].level);
+ block->getKernel()->setDataSet(dataSetPtr);
+ }
+
+ delete [] dataSetArray;
+}
+
+void MPIIORestartCoProcessor::readBoundaryConds()
+{
+ int rank, size;
+ MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+ MPI_Comm_size(MPI_COMM_WORLD, &size);
+
+ MPI_File file_handler;
+ int error = MPI_File_open(MPI_COMM_WORLD, "outputBoundCond.bin", MPI_MODE_CREATE | MPI_MODE_RDWR , MPI_INFO_NULL, &file_handler);
+
+ int blocksCount = 0;
+ int dataCount1000 = 0;
+ int dataCount2 = 0;
+ // read count of blocks
+ MPI_File_read_at(file_handler, rank * sizeof(int), &blocksCount, 1, MPI_INT, MPI_STATUS_IGNORE);
+ // read count of big BoundaryCondition blocks
+ MPI_File_read_at(file_handler, (rank + size) * sizeof(int), &dataCount1000, 1, MPI_INT, MPI_STATUS_IGNORE);
+ // read count of indexContainer values in all blocks
+ MPI_File_read_at(file_handler, (rank + 2 * size) * sizeof(int), &dataCount2, 1, MPI_INT, MPI_STATUS_IGNORE);
+
+ size_t dataCount = dataCount1000 * BLOCK_SIZE;
+ BCAdd* bcAddArray = new BCAdd[blocksCount];
+ BoundaryCondition* bcArray = new BoundaryCondition[dataCount];
+ BoundaryCondition* nullBouCond = new BoundaryCondition();
+ memset(nullBouCond, 0, sizeof(BoundaryCondition));
+ int* intArray1 = new int[blocksCount * blockParamStr.bcindexmatrix_count];
+ int* intArray2 = new int[dataCount2];
+
+ size_t read_offset = 3 * size * sizeof(int);
+ size_t next_read_offset = 0;
+
+ if(size > 1)
+ {
+ if(rank == 0)
+ {
+ next_read_offset = read_offset + blocksCount * sizeof(BCAdd) + dataCount * sizeof(BoundaryCondition) + (blocksCount * blockParamStr.bcindexmatrix_count + dataCount2) * sizeof(int);
+ MPI_Send(&next_read_offset, 1, MPI_LONG_LONG_INT, 1, 5, MPI_COMM_WORLD);
+ }
+ else
+ {
+ MPI_Recv(&read_offset, 1, MPI_LONG_LONG_INT, rank - 1, 5, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
+ next_read_offset = read_offset + blocksCount * sizeof(BCAdd) + dataCount * sizeof(BoundaryCondition) + (blocksCount * blockParamStr.bcindexmatrix_count + dataCount2) * sizeof(int);
+ if(rank < size - 1)
+ MPI_Send(&next_read_offset, 1, MPI_LONG_LONG_INT, rank + 1, 5, MPI_COMM_WORLD);
+ }
+ }
+
+ MPI_File_read_at(file_handler, read_offset, &bcAddArray[0], blocksCount, boundCondTypeAdd, MPI_STATUS_IGNORE);
+ MPI_File_read_at(file_handler, read_offset + blocksCount * sizeof(BCAdd), &bcArray[0], dataCount1000, boundCondType1000, MPI_STATUS_IGNORE);
+ MPI_File_read_at(file_handler, read_offset + blocksCount * sizeof(BCAdd) + dataCount * sizeof(BoundaryCondition), &intArray1[0], blocksCount, bcindexmatrixType, MPI_STATUS_IGNORE);
+ MPI_File_read_at(file_handler, read_offset + blocksCount * sizeof(BCAdd) + dataCount * sizeof(BoundaryCondition) + blocksCount * blockParamStr.bcindexmatrix_count * sizeof(int), &intArray2[0], dataCount2, MPI_INT, MPI_STATUS_IGNORE);
+ MPI_File_sync(file_handler);
+
+ MPI_File_close(&file_handler);
+
+ int index = 0, index1 = 0, index2 = 0;
+ std::vector<BoundaryConditionsPtr> bcVector;
+ std::vector<int> bcindexmatrixV;
+ std::vector<int> indexContainerV;
+
+ for(size_t n = 0; n < blocksCount; n++)
+ {
+ bcVector.resize(0);
+ bcindexmatrixV.resize(0);
+ indexContainerV.resize(0);
+
+ for(size_t ibc = 0; ibc < bcAddArray[n].boundCond_count; ibc++)
+ {
+ BoundaryConditionsPtr bc;
+ if(memcmp(&bcArray[index], nullBouCond, sizeof(BoundaryCondition)) == 0)
+ bc = BoundaryConditionsPtr();
+ else
+ {
+ bc = BoundaryConditionsPtr(new BoundaryConditions);
+ bc->noslipBoundaryFlags = bcArray[index].noslipBoundaryFlags;
+ bc->slipBoundaryFlags = bcArray[index].slipBoundaryFlags;
+ bc->densityBoundaryFlags = bcArray[index].densityBoundaryFlags;
+ bc->velocityBoundaryFlags = bcArray[index].velocityBoundaryFlags;
+ bc->wallModelBoundaryFlags = bcArray[index].wallModelBoundaryFlags;
+ bc->bcVelocityX1 = bcArray[index].bcVelocityX1;
+ bc->bcVelocityX2 = bcArray[index].bcVelocityX2;
+ bc->bcVelocityX3 = bcArray[index].bcVelocityX3;
+ bc->bcDensity = bcArray[index].bcDensity;
+ bc->bcLodiDensity = bcArray[index].bcLodiDensity;
+ bc->bcLodiVelocityX1 = bcArray[index].bcLodiVelocityX1;
+ bc->bcLodiVelocityX2 = bcArray[index].bcLodiVelocityX2;
+ bc->bcLodiVelocityX3 = bcArray[index].bcLodiVelocityX3;
+ bc->bcLodiLentgh = bcArray[index].bcLodiLentgh;
+
+ bc->nx1 = bcArray[index].nx1;
+ bc->nx2 = bcArray[index].nx2;
+ bc->nx3 = bcArray[index].nx3;
+ for(int iq = 0; iq < 26; iq++)
+ bc->setQ(bcArray[index].q[iq], iq);
+ bc->setBcAlgorithmType(bcArray[index].algorithmType);
+ }
+
+ bcVector.push_back(bc);
+ index++;
+ }
+
+ for(int b1 = 0; b1 < blockParamStr.bcindexmatrix_count; b1++)
+ bcindexmatrixV.push_back(intArray1[index1++]);
+
+ for(int b2 = 0; b2 < bcAddArray[n].indexContainer_count; b2++)
+ indexContainerV.push_back(intArray2[index2++]);
+
+ CbArray3D<int,IndexerX3X2X1> bcim(bcindexmatrixV, blockParamStr.nx1, blockParamStr.nx2, blockParamStr.nx3);
+
+ BCProcessorPtr bcProc(new BCProcessor());
+ BCArray3D &bcArr = bcProc->getBCArray();
+ bcArr.bcindexmatrix = bcim;
+ bcArr.bcvector = bcVector;
+ bcArr.indexContainer = indexContainerV;
+
+ Block3DPtr block = grid->getBlock(bcAddArray[n].x1, bcAddArray[n].x2, bcAddArray[n].x3, bcAddArray[n].level);
+ block->getKernel()->setBCProcessor(bcProc);
+ }
+
+ delete nullBouCond;
+ delete[] bcArray;
+ delete[] bcAddArray;
+ delete[] intArray1;
+ delete[] intArray2;
+}
+
diff --git a/source/VirtualFluidsCore/CoProcessors/MPIIORestartCoProcessor.h b/source/VirtualFluidsCore/CoProcessors/MPIIORestartCoProcessor.h
new file mode 100644
index 0000000000000000000000000000000000000000..a626c364205e04c46ebb7cb31c02eeccae074786
--- /dev/null
+++ b/source/VirtualFluidsCore/CoProcessors/MPIIORestartCoProcessor.h
@@ -0,0 +1,164 @@
+#ifndef _MPIWRITEBLOCKSCOPROCESSOR_H_
+#define _MPIWRITEBLOCKSCOPROCESSOR_H_
+
+#include "mpi.h"
+
+#include "CoProcessor.h"
+#include "Communicator.h"
+#include "WbWriter.h"
+
+#include <boost/shared_ptr.hpp>
+
+class MPIIORestartCoProcessor;
+typedef boost::shared_ptr<MPIIORestartCoProcessor> MPIIORestartCoProcessorPtr;
+
+//! \class MPIWriteBlocksCoProcessor
+//! \brief Writes the grid each timestep into the files and reads the grip from the files before regenerating
+class MPIIORestartCoProcessor: public CoProcessor
+{
+ //! \struct GridParam
+ //! \brief Structure describes parameters of the grid
+ //! \details The structure is nessasary to restore the grid correctly
+ typedef struct GridParam
+ {
+ double trafoParams[33];
+ double deltaX;
+ int blockNx1;
+ int blockNx2;
+ int blockNx3;
+ int nx1;
+ int nx2;
+ int nx3;
+ bool periodicX1;
+ bool periodicX2;
+ bool periodicX3;
+ bool active;
+ bool transformation;
+ }GridParam;
+
+ //! \struct blockParam
+ //! \brief Structure describes parameters of the block that are equal in all blocks
+ //! \details The structure used to store some parameters needed to restore dataSet arrays
+ typedef struct blockParam
+ {
+ int nx1; // to find the right block
+ int nx2;
+ int nx3;
+ int nx[9][4]; // 9 arrays x (nx1, nx2, nx3, nx4)
+ int doubleCountInBlock; // how many double-values are in all arrays dataSet in one (any) block
+ int bcindexmatrix_count; // how many bcindexmatrix-values are in one (any) block
+ }blockParam;
+
+ //! \struct Block3d
+ //! \brief Structure contains information of the block
+ //! \details The structure is used to write the data describing the block in the grid when saving the grid
+ //! and to read it when restoring the grid
+ typedef struct Block3d
+ {
+ double collFactor;
+ double deltaT;
+ int x1;
+ int x2;
+ int x3;
+ int bundle;
+ int rank;
+ int lrank;
+ int part;
+ int globalID;
+ int localID;
+ int level;
+ int interpolationFlagCF;
+ int interpolationFlagFC;
+ int counter;
+ int ghostLayerWidth;
+ bool active;
+ bool compressible;
+ bool withForcing;
+ }Block3d;
+
+ //! \struct dataSet
+ //! \brief Structure containes information identifying the block
+ //! \details The structure is used to find the needed block in the grid when restoring a dataSet
+ typedef struct dataSet
+ {
+ int x1;
+ int x2;
+ int x3;
+ int level;
+ }dataSet;
+
+ //! \struct BoundaryCondition
+ //! \brief Structure containes information about boundary conditions of the block
+ //! \details The structure is used to write data describing boundary conditions of the blocks when saving the grid
+ //! and to read it when restoring the grid
+ typedef struct BoundaryCondition
+ {
+ long long noslipBoundaryFlags; // MPI_LONG_LONG
+ long long slipBoundaryFlags;
+ long long velocityBoundaryFlags;
+ long long densityBoundaryFlags;
+ long long wallModelBoundaryFlags;
+
+ float bcVelocityX1;
+ float bcVelocityX2;
+ float bcVelocityX3;
+ float bcDensity;
+
+ float bcLodiDensity;
+ float bcLodiVelocityX1;
+ float bcLodiVelocityX2;
+ float bcLodiVelocityX3;
+ float bcLodiLentgh;
+
+ float nx1,nx2,nx3;
+ float q[26]; // MPI_FLOAT
+
+ char algorithmType;
+ }BoundaryCondition;
+
+ //! \struct BCAdd
+ //! \brief Structure containes information identifying the block
+ //! and some parameters of the arrays of boundary conditions that are equal in all blocks
+ //! \details The structure is used to find the needed block in the grid when restoring a dataSet
+ //! and to set common parameters
+ typedef struct BCAdd
+ {
+ int x1; // to find the right block
+ int x2;
+ int x3;
+ int level;
+ int boundCond_count; // how many BoundaryCondition-structures are in this block
+ int indexContainer_count; // how many indexContainer-values are in this block
+ }BCAdd;
+
+public:
+ MPIIORestartCoProcessor(Grid3DPtr grid, UbSchedulerPtr s, const std::string& path, CommunicatorPtr comm);
+ virtual ~MPIIORestartCoProcessor();
+ //! Each timestep writes the grid into the files
+ void process(double step);
+ //! Reads the grid from the files before grid reconstruction
+ void restart();
+ //! Writes the blocks of the grid into the file outputBlocks.bin
+ void writeBlocks();
+ //! Writes the datasets of the blocks into the file outputDataSet.bin
+ void writeDataSet();
+ //! Writes the boundary conditions of the blocks into the file outputBoundCond.bin
+ void writeBoundaryConds();
+ //! Reads the blocks of the grid from the file outputBlocks.bin
+ void readBlocks();
+ //! Reads the datasets of the blocks from the file outputDataSet.bin
+ void readDataSet();
+ //! Reads the boundary conditions of the blocks from the file outputBoundCond.bin
+ void readBoundaryConds();
+
+protected:
+ std::string path;
+ CommunicatorPtr comm;
+
+private:
+ MPI_Datatype gridParamType, blockParamType, block3dType, dataSetType, dataSetDoubleType, boundCondType, boundCondType1000, boundCondTypeAdd, bcindexmatrixType;
+ blockParam blockParamStr;
+
+};
+
+#endif
diff --git a/source/VirtualFluidsCore/Data/D3Q27EsoTwist3DSplittedVector.cpp b/source/VirtualFluidsCore/Data/D3Q27EsoTwist3DSplittedVector.cpp
index e8a1185fa59f1d98f9267b39d6cedb10d2562f96..d048f9c58ab8c08492807ed6b16c5051ad8b5145 100644
--- a/source/VirtualFluidsCore/Data/D3Q27EsoTwist3DSplittedVector.cpp
+++ b/source/VirtualFluidsCore/Data/D3Q27EsoTwist3DSplittedVector.cpp
@@ -619,9 +619,34 @@ CbArray3D<LBMReal,IndexerX3X2X1>::CbArray3DPtr D3Q27EsoTwist3DSplittedVector::ge
return this->zeroDistributions;
}
//////////////////////////////////////////////////////////////////////////
-void D3Q27EsoTwist3DSplittedVector::getDistributionAfterLastStep( LBMReal* const f, size_t x1, size_t x2, size_t x3 )
+void D3Q27EsoTwist3DSplittedVector::setNX1(size_t newNX1)
{
-
+ NX1 = newNX1;
+}
+//////////////////////////////////////////////////////////////////////////
+void D3Q27EsoTwist3DSplittedVector::setNX2(size_t newNX2)
+{
+ NX2 = newNX2;
+}
+//////////////////////////////////////////////////////////////////////////
+void D3Q27EsoTwist3DSplittedVector::setNX3(size_t newNX3)
+{
+ NX3 = newNX3;
+}
+//////////////////////////////////////////////////////////////////////////
+void D3Q27EsoTwist3DSplittedVector::setLocalDistributions(CbArray4D<LBMReal, IndexerX4X3X2X1>::CbArray4DPtr array)
+{
+ localDistributions = array;
+}
+//////////////////////////////////////////////////////////////////////////
+void D3Q27EsoTwist3DSplittedVector::setNonLocalDistributions(CbArray4D<LBMReal, IndexerX4X3X2X1>::CbArray4DPtr array)
+{
+ nonLocalDistributions = array;
+}
+//////////////////////////////////////////////////////////////////////////
+void D3Q27EsoTwist3DSplittedVector::setZeroDistributions(CbArray3D<LBMReal, IndexerX3X2X1>::CbArray3DPtr array)
+{
+ zeroDistributions = array;
}
//////////////////////////////////////////////////////////////////////////
diff --git a/source/VirtualFluidsCore/Data/D3Q27EsoTwist3DSplittedVector.h b/source/VirtualFluidsCore/Data/D3Q27EsoTwist3DSplittedVector.h
index 126c0f54a33237dc48c16e785eb1a070b84e0667..c3badacb856af044035544409ba5e477a708b9ca 100644
--- a/source/VirtualFluidsCore/Data/D3Q27EsoTwist3DSplittedVector.h
+++ b/source/VirtualFluidsCore/Data/D3Q27EsoTwist3DSplittedVector.h
@@ -50,7 +50,12 @@ public:
//////////////////////////////////////////////////////////////////////////
CbArray3D<LBMReal,IndexerX3X2X1>::CbArray3DPtr getZeroDistributions();
//////////////////////////////////////////////////////////////////////////
- void getDistributionAfterLastStep(LBMReal* const f, size_t x1, size_t x2, size_t x3);
+ void setNX1(size_t newNX1);
+ void setNX2(size_t newNX2);
+ void setNX3(size_t newNX3);
+ void setLocalDistributions(CbArray4D<LBMReal, IndexerX4X3X2X1>::CbArray4DPtr array);
+ void setNonLocalDistributions(CbArray4D<LBMReal, IndexerX4X3X2X1>::CbArray4DPtr array);
+ void setZeroDistributions(CbArray3D<LBMReal, IndexerX3X2X1>::CbArray3DPtr array);
protected:
CbArray4D<LBMReal,IndexerX4X3X2X1>::CbArray4DPtr localDistributions;
@@ -58,6 +63,8 @@ protected:
CbArray3D<LBMReal,IndexerX3X2X1>::CbArray3DPtr zeroDistributions;
size_t NX1, NX2, NX3;
+ friend class MPIIORestartCoProcessor;
+
friend class boost::serialization::access;
template<class Archive>
void serialize(Archive & ar, const unsigned int version)
diff --git a/source/VirtualFluidsCore/Grid/Block3D.h b/source/VirtualFluidsCore/Grid/Block3D.h
index dc9cbb1038ad192890a22c73204601b324f74215..c0ce971c0dbb67f2cd25b760577cb4b3db9ca4eb 100644
--- a/source/VirtualFluidsCore/Grid/Block3D.h
+++ b/source/VirtualFluidsCore/Grid/Block3D.h
@@ -131,6 +131,8 @@ private:
int level;
static int counter;
+ friend class MPIIORestartCoProcessor;
+
friend class boost::serialization::access;
template<class Archive>
void serialize(Archive & ar, const unsigned int version)
diff --git a/source/VirtualFluidsCore/LBM/LBMKernel.cpp b/source/VirtualFluidsCore/LBM/LBMKernel.cpp
index d301691eef7278915e62bd13822bac4edf6b0acd..b3636973238d9d5777e0c0a210a8ecb99a992a5d 100644
--- a/source/VirtualFluidsCore/LBM/LBMKernel.cpp
+++ b/source/VirtualFluidsCore/LBM/LBMKernel.cpp
@@ -133,6 +133,11 @@ DataSet3DPtr LBMKernel::getDataSet() const
return this->dataSet;
}
//////////////////////////////////////////////////////////////////////////
+LBMReal LBMKernel::getDeltaT()
+{
+ return this->deltaT;
+}
+//////////////////////////////////////////////////////////////////////////
void LBMKernel::setDeltaT( LBMReal dt )
{
deltaT = dt;
diff --git a/source/VirtualFluidsCore/LBM/LBMKernel.h b/source/VirtualFluidsCore/LBM/LBMKernel.h
index c36a5229ca84963f7ec8b1ec3aa40af23ee37006..c88a7f92a9a2af285031a314eb07a2d3f9cd2565 100644
--- a/source/VirtualFluidsCore/LBM/LBMKernel.h
+++ b/source/VirtualFluidsCore/LBM/LBMKernel.h
@@ -58,6 +58,7 @@ public:
void setIndex(int x1, int x2, int x3);
+ LBMReal getDeltaT();
void setDeltaT(LBMReal dt);
bool getCompressible() const;