diff --git a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C index 8a2bb3ed8..1bb72b7f5 100644 --- a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C +++ b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C @@ -119,6 +119,7 @@ int main(int argc, char *argv[]) bool writeCellDist = args.optionFound("cellDist"); bool copyUniform = args.optionFound("copyUniform"); + bool decomposeMeshOnly = args.optionFound("mesh"); bool decomposeFieldsOnly = args.optionFound("fields"); bool filterPatches = args.optionFound("filterPatches"); bool forceOverwrite = args.optionFound("force"); @@ -136,9 +137,9 @@ int main(int argc, char *argv[]) ( runTime.path() /(word("processor") + name(nProcs)) - /runTime.constant() - /regionDir - /polyMesh::meshSubDir + ///runTime.constant() + ///regionDir + ///polyMesh::meshSubDir ) ) { @@ -263,7 +264,6 @@ int main(int argc, char *argv[]) ) ); - // Decompose the mesh if (!decomposeFieldsOnly) { @@ -271,61 +271,11 @@ int main(int argc, char *argv[]) meshDecomp.writeDecomposition(); - if (writeCellDist) - { - const labelList& procIds = meshDecomp.cellToProc(); - - // Write the decomposition as labelList for use with 'manual' - // decomposition method. - labelIOList cellDecomposition - ( - IOobject - ( - "cellDecomposition", - mesh.facesInstance(), - mesh, - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - procIds - ); - cellDecomposition.write(); - - Info<< nl << "Wrote decomposition to " - << cellDecomposition.objectPath() - << " for use in manual decomposition." << endl; - - // Write as volScalarField for post-processing - volScalarField cellDist - ( - IOobject - ( - "cellDist", - runTime.timeName(), - mesh.dbDir(), - mesh, - IOobject::NO_READ, - IOobject::NO_WRITE - ), - mesh, - dimensionedScalar("cellDist", dimless, 0), - zeroGradientFvPatchScalarField::typeName - ); - - forAll(procIds, celli) - { - cellDist[celli] = procIds[celli]; - } - - cellDist.write(); - - Info<< nl << "Wrote decomposition as volScalarField to " - << cellDist.name() << " for use in post-processing." - << endl; - } + SliceStreamRepo::instance()->clear(); } + if (!decomposeMeshOnly) + { // Search for list of objects for this time IOobjectList objects(mesh, runTime.timeName()); @@ -361,242 +311,6 @@ int main(int argc, char *argv[]) PtrList surfaceTensorFields; readFields(mesh, objects, surfaceTensorFields); - - // Construct the point fields - // ~~~~~~~~~~~~~~~~~~~~~~~~~~ - const pointMesh& pMesh = pointMesh::New(mesh); - - PtrList pointScalarFields; - readFields(pMesh, objects, pointScalarFields); - - PtrList pointVectorFields; - readFields(pMesh, objects, pointVectorFields); - - PtrList pointSphericalTensorFields; - readFields(pMesh, objects, pointSphericalTensorFields); - - PtrList pointSymmTensorFields; - readFields(pMesh, objects, pointSymmTensorFields); - - PtrList pointTensorFields; - readFields(pMesh, objects, pointTensorFields); - - - // Construct the tetPoint fields - // ~~~~~~~~~~~~~~~~~~~~~~~~~~ - tetPolyMesh* tetMeshPtr = nullptr; - - PtrList tetPointScalarFields; - PtrList tetPointVectorFields; - PtrList tetPointSphericalTensorFields; - PtrList tetPointSymmTensorFields; - PtrList tetPointTensorFields; - - PtrList elementScalarFields; - PtrList elementVectorFields; - - if - ( - objects.lookupClass("tetPointScalarField").size() > 0 - || objects.lookupClass("tetPointVectorField").size() > 0 - || objects.lookupClass("tetPointSphericalTensorField").size() > 0 - || objects.lookupClass("tetPointSymmTensorField").size() > 0 - || objects.lookupClass("tetPointTensorField").size() > 0 - - || objects.lookupClass("elementScalarField").size() > 0 - || objects.lookupClass("elementVectorField").size() > 0 - ) - { - tetMeshPtr = new tetPolyMesh(mesh); - tetPolyMesh& tetMesh = *tetMeshPtr; - - readFields(tetMesh, objects, tetPointScalarFields); - readFields(tetMesh, objects, tetPointVectorFields); - readFields(tetMesh, objects, tetPointSphericalTensorFields); - readFields(tetMesh, objects, tetPointSymmTensorFields); - readFields(tetMesh, objects, tetPointTensorFields); - - readFields(tetMesh, objects, elementScalarFields); - readFields(tetMesh, objects, elementVectorFields); - } - - - // Construct the Lagrangian fields - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - - fileNameList cloudDirs - ( - readDir(runTime.timePath()/cloud::prefix, fileName::DIRECTORY) - ); - - // Particles - PtrList > lagrangianPositions(cloudDirs.size()); - // Particles per cell - PtrList< List*> > cellParticles(cloudDirs.size()); - - PtrList > lagrangianLabelFields(cloudDirs.size()); - PtrList > lagrangianScalarFields(cloudDirs.size()); - PtrList > lagrangianVectorFields(cloudDirs.size()); - PtrList > lagrangianSphericalTensorFields - ( - cloudDirs.size() - ); - PtrList > lagrangianSymmTensorFields - ( - cloudDirs.size() - ); - PtrList > lagrangianTensorFields - ( - cloudDirs.size() - ); - - label cloudI = 0; - - forAll(cloudDirs, i) - { - IOobjectList sprayObjs - ( - mesh, - runTime.timeName(), - cloud::prefix/cloudDirs[i] - ); - - IOobject* positionsPtr = sprayObjs.lookup("positions"); - - if (positionsPtr) - { - // Read lagrangian particles - // ~~~~~~~~~~~~~~~~~~~~~~~~~ - - Info<< "Identified lagrangian data set: " << cloudDirs[i] << endl; - - lagrangianPositions.set - ( - cloudI, - new Cloud - ( - mesh, - cloudDirs[i], - false - ) - ); - - - // Sort particles per cell - // ~~~~~~~~~~~~~~~~~~~~~~~ - - cellParticles.set - ( - cloudI, - new List*> - ( - mesh.nCells(), - static_cast*>(nullptr) - ) - ); - - label i = 0; - - forAllIter - ( - Cloud, - lagrangianPositions[cloudI], - iter - ) - { - iter().index() = i++; - - label celli = iter().cell(); - - // Check - if (celli < 0 || celli >= mesh.nCells()) - { - FatalErrorIn(args.executable()) - << "Illegal cell number " << celli - << " for particle with index " << iter().index() - << " at position " << iter().position() << nl - << "Cell number should be between 0 and " - << mesh.nCells()-1 << nl - << "On this mesh the particle should be in cell " - << mesh.findCell(iter().position()) - << exit(FatalError); - } - - if (!cellParticles[cloudI][celli]) - { - cellParticles[cloudI][celli] = - new SLList(); - } - - cellParticles[cloudI][celli]->append(&iter()); - } - - // Read fields - // ~~~~~~~~~~~ - - IOobjectList lagrangianObjects - ( - mesh, - runTime.timeName(), - cloud::prefix/cloudDirs[cloudI] - ); - - lagrangianFieldDecomposer::readFields - ( - cloudI, - lagrangianObjects, - lagrangianLabelFields - ); - - lagrangianFieldDecomposer::readFields - ( - cloudI, - lagrangianObjects, - lagrangianScalarFields - ); - - lagrangianFieldDecomposer::readFields - ( - cloudI, - lagrangianObjects, - lagrangianVectorFields - ); - - lagrangianFieldDecomposer::readFields - ( - cloudI, - lagrangianObjects, - lagrangianSphericalTensorFields - ); - - lagrangianFieldDecomposer::readFields - ( - cloudI, - lagrangianObjects, - lagrangianSymmTensorFields - ); - - lagrangianFieldDecomposer::readFields - ( - cloudI, - lagrangianObjects, - lagrangianTensorFields - ); - - cloudI++; - } - } - - lagrangianPositions.setSize(cloudI); - cellParticles.setSize(cloudI); - lagrangianLabelFields.setSize(cloudI); - lagrangianScalarFields.setSize(cloudI); - lagrangianVectorFields.setSize(cloudI); - lagrangianSphericalTensorFields.setSize(cloudI); - lagrangianSymmTensorFields.setSize(cloudI); - lagrangianTensorFields.setSize(cloudI); - - // Any uniform data to copy/link? fileName uniformDir("uniform"); @@ -613,6 +327,72 @@ int main(int argc, char *argv[]) Info<< endl; + if ( mesh.time().writeFormat() == IOstream::COHERENT ) { + + Time myDb( Time::controlDictName, + mesh.time().rootPath(), + mesh.time().caseName(), + "system", + "constant", + true ); + + // Read the mesh + fvMesh coherentMesh + ( + IOobject + ( + regionName, + myDb.timeName(), + myDb + ) + ); + + auto coherentFaceAddressing = meshDecomp.coherentFaceAddressing(); + auto coherentCellAddressing = meshDecomp.coherentCellAddressing(); + auto coherentBoundaryAddressing = meshDecomp.blankAddressing( mesh.boundaryMesh().size() ); + + // Old format of face addressing requires a plus one. + for ( auto& faceId : coherentFaceAddressing ) { ++faceId; } + + // FV fields + if + ( + volScalarFields.size() + || volVectorFields.size() + || volSphericalTensorFields.size() + || volSymmTensorFields.size() + || volTensorFields.size() + || surfaceScalarFields.size() + || surfaceVectorFields.size() + || surfaceSphericalTensorFields.size() + || surfaceSymmTensorFields.size() + || surfaceTensorFields.size() + ) + { + fvFieldDecomposer fieldDecomposer + ( + mesh, + coherentMesh, + coherentFaceAddressing, + coherentCellAddressing, + coherentBoundaryAddressing + ); + + fieldDecomposer.decomposeFields(volScalarFields); + fieldDecomposer.decomposeFields(volVectorFields); + fieldDecomposer.decomposeFields(volSphericalTensorFields); + fieldDecomposer.decomposeFields(volSymmTensorFields); + fieldDecomposer.decomposeFields(volTensorFields); + + fieldDecomposer.decomposeFields(surfaceScalarFields); + fieldDecomposer.decomposeFields(surfaceVectorFields); + fieldDecomposer.decomposeFields(surfaceSphericalTensorFields); + fieldDecomposer.decomposeFields(surfaceSymmTensorFields); + fieldDecomposer.decomposeFields(surfaceTensorFields); + } + } + else + { // Split the fields over processors for (label procI = 0; procI < meshDecomp.nProcs(); procI++) { @@ -741,123 +521,6 @@ int main(int argc, char *argv[]) fieldDecomposer.decomposeFields(surfaceTensorFields); } - - // Point fields - if - ( - pointScalarFields.size() - || pointVectorFields.size() - || pointSphericalTensorFields.size() - || pointSymmTensorFields.size() - || pointTensorFields.size() - ) - { - const pointMesh& procPMesh = pointMesh::New(procMesh, true); - - pointFieldDecomposer fieldDecomposer - ( - pMesh, - procPMesh, - pointProcAddressing, - boundaryProcAddressing - ); - - fieldDecomposer.decomposeFields(pointScalarFields); - fieldDecomposer.decomposeFields(pointVectorFields); - fieldDecomposer.decomposeFields(pointSphericalTensorFields); - fieldDecomposer.decomposeFields(pointSymmTensorFields); - fieldDecomposer.decomposeFields(pointTensorFields); - } - - - // tetPoint fields - if (tetMeshPtr) - { - const tetPolyMesh& tetMesh = *tetMeshPtr; - tetPolyMesh procTetMesh(procMesh); - - tetPointFieldDecomposer fieldDecomposer - ( - tetMesh, - procTetMesh, - pointProcAddressing, - faceProcAddressing, - cellProcAddressing, - boundaryProcAddressing - ); - - fieldDecomposer.decomposeFields(tetPointScalarFields); - fieldDecomposer.decomposeFields(tetPointVectorFields); - fieldDecomposer.decomposeFields(tetPointSphericalTensorFields); - fieldDecomposer.decomposeFields(tetPointSymmTensorFields); - fieldDecomposer.decomposeFields(tetPointTensorFields); - - fieldDecomposer.decomposeFields(elementScalarFields); - fieldDecomposer.decomposeFields(elementVectorFields); - } - - - // If there is lagrangian data write it out - forAll(lagrangianPositions, cloudI) - { - if (lagrangianPositions[cloudI].size()) - { - lagrangianFieldDecomposer fieldDecomposer - ( - mesh, - procMesh, - cellProcAddressing, - cloudDirs[cloudI], - lagrangianPositions[cloudI], - cellParticles[cloudI] - ); - - // Lagrangian fields - if - ( - lagrangianLabelFields[cloudI].size() - || lagrangianScalarFields[cloudI].size() - || lagrangianVectorFields[cloudI].size() - || lagrangianSphericalTensorFields[cloudI].size() - || lagrangianSymmTensorFields[cloudI].size() - || lagrangianTensorFields[cloudI].size() - ) - { - fieldDecomposer.decomposeFields - ( - cloudDirs[cloudI], - lagrangianLabelFields[cloudI] - ); - fieldDecomposer.decomposeFields - ( - cloudDirs[cloudI], - lagrangianScalarFields[cloudI] - ); - fieldDecomposer.decomposeFields - ( - cloudDirs[cloudI], - lagrangianVectorFields[cloudI] - ); - fieldDecomposer.decomposeFields - ( - cloudDirs[cloudI], - lagrangianSphericalTensorFields[cloudI] - ); - fieldDecomposer.decomposeFields - ( - cloudDirs[cloudI], - lagrangianSymmTensorFields[cloudI] - ); - fieldDecomposer.decomposeFields - ( - cloudDirs[cloudI], - lagrangianTensorFields[cloudI] - ); - } - } - } - - // Any non-decomposed data to copy? if (uniformDir.size()) { @@ -889,176 +552,10 @@ int main(int argc, char *argv[]) } } } + } - - - if (tetMeshPtr) - { - delete tetMeshPtr; - tetMeshPtr = nullptr; - } - - - // Finite area mesh and field decomposition - - IOobject faMeshBoundaryIOobj - ( - "faBoundary", - mesh.time().findInstance - ( - mesh.meshDir(), - "boundary" - ), - faMesh::meshSubDir, - mesh, - IOobject::READ_IF_PRESENT, - IOobject::NO_WRITE - ); - - - if(faMeshBoundaryIOobj.headerOk()) - { - Info << "\nFinite area mesh decomposition" << endl; - - faMeshDecomposition aMesh(mesh); - - aMesh.decomposeMesh(filterPatches); - - aMesh.writeDecomposition(); - - - // Construct the area fields - // ~~~~~~~~~~~~~~~~~~~~~~~~ - PtrList areaScalarFields; - readFields(aMesh, objects, areaScalarFields); - - PtrList areaVectorFields; - readFields(aMesh, objects, areaVectorFields); - - PtrList areaSphericalTensorFields; - readFields(aMesh, objects, areaSphericalTensorFields); - - PtrList areaSymmTensorFields; - readFields(aMesh, objects, areaSymmTensorFields); - - PtrList areaTensorFields; - readFields(aMesh, objects, areaTensorFields); - - - // Construct the edge fields - // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - PtrList edgeScalarFields; - readFields(aMesh, objects, edgeScalarFields); - - Info << endl; - - // Split the fields over processors - for (label procI = 0; procI < meshDecomp.nProcs(); procI++) - { - Info<< "Processor " << procI - << ": finite area field transfer" << endl; - - // open the database - Time processorDb - ( - Time::controlDictName, - args.rootPath(), - args.caseName()/fileName(word("processor") + name(procI)) - ); - - processorDb.setTime(runTime); - - // Read the mesh - fvMesh procFvMesh - ( - IOobject - ( - regionName, - processorDb.timeName(), - processorDb - ) - ); - - faMesh procMesh(procFvMesh); - - word faceLabelsInstance = - procMesh.time().findInstance - ( - procMesh.meshDir(), - "faceLabels" - ); - - - labelIOList faceProcAddressing - ( - IOobject - ( - "faceProcAddressing", - faceLabelsInstance, - procMesh.meshSubDir, - procFvMesh, - IOobject::MUST_READ, - IOobject::NO_WRITE - ) - ); - - labelIOList boundaryProcAddressing - ( - IOobject - ( - "boundaryProcAddressing", - faceLabelsInstance, - procMesh.meshSubDir, - procFvMesh, - IOobject::MUST_READ, - IOobject::NO_WRITE - ) - ); - - // FA fields - if - ( - areaScalarFields.size() - || areaVectorFields.size() - || areaSphericalTensorFields.size() - || areaSymmTensorFields.size() - || areaTensorFields.size() - || edgeScalarFields.size() - ) - { - labelIOList edgeProcAddressing - ( - IOobject - ( - "edgeProcAddressing", - faceLabelsInstance, - procMesh.meshSubDir, - procFvMesh, - IOobject::MUST_READ, - IOobject::NO_WRITE - ) - ); - - faFieldDecomposer fieldDecomposer - ( - aMesh, - procMesh, - edgeProcAddressing, - faceProcAddressing, - boundaryProcAddressing - ); - - fieldDecomposer.decomposeFields(areaScalarFields); - fieldDecomposer.decomposeFields(areaVectorFields); - fieldDecomposer.decomposeFields(areaSphericalTensorFields); - fieldDecomposer.decomposeFields(areaSymmTensorFields); - fieldDecomposer.decomposeFields(areaTensorFields); - - fieldDecomposer.decomposeFields(edgeScalarFields); - } - } - } - + } // for-loop on procI + } // if ( IOstream::COHERENT ) Info<< "\nEnd.\n" << endl; diff --git a/etc/bashrc b/etc/bashrc index 9ddc4dcac..cf95a07cb 100755 --- a/etc/bashrc +++ b/etc/bashrc @@ -34,36 +34,35 @@ export WM_PROJECT=foam export WM_FORK=extend -export WM_PROJECT_VERSION=4.1 +export WM_PROJECT_VERSION=4.1r # helps to easily write #ifdefs to detect a dev-version export FOAM_DEV=1 -############################################################################### -# USER EDITABLE PART +#------------------------------------------------------------------------------ +# (advanced / legacy) + +# [projectDir] - directory containing this OpenFOAM version. +# \- When this file is located as $WM_PROJECT_DIR/etc/bashrc, the next lines +# should work when sourced by BASH or ZSH shells. If this however fails, +# set one of the fallback values to an appropriate path. # -# either set $FOAM_INST_DIR before sourcing this file or set -# $foamInstall below to where FOAM is installed -# -# Location of FOAM installation -# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -foamInstall=$HOME/$WM_PROJECT -# foamInstall=~$WM_PROJECT -# foamInstall=/usr/local/$WM_PROJECT -# foamInstall=/opt/$WM_PROJECT -# foamInstall=/usr/lib -# -# END OF (NORMAL) USER EDITABLE PART +# This can be removed if an absolute path is provided for WM_PROJECT_DIR +# later on in this file +# -- +projectDir="${BASH_SOURCE:-${ZSH_NAME:+$0}}"; +[ -n "$projectDir" ] && projectDir="$(\cd $(dirname $projectDir)/.. && \pwd -L)" ||\ +projectDir="$HOME/openfoam/foam-extend-$WM_PROJECT_VERSION" ################################################################################ # Set $USER if it has not been set (eg. on Windows/MSYS systems) -if [ -z $USER ] +if [ -z "$USER" ] then - export USER=`whoami` + export USER="$(whoami)" fi # note the location for later use (eg, in job scripts) -: ${FOAM_INST_DIR:=$foamInstall}; export FOAM_INST_DIR +: ${FOAM_INST_DIR:=$(dirname $projectDir)}; export FOAM_INST_DIR # The old dirs to be cleaned from the various environment variables # - remove anything under top-level directory. @@ -74,16 +73,17 @@ then foamOldDirs="$WM_PROJECT_INST_DIR $foamOldDirs" fi - # Location of site/user files # ~~~~~~~~~~~~~~~~~~~~~~~~~~~ export WM_PROJECT_INST_DIR=$FOAM_INST_DIR -export WM_PROJECT_DIR=$WM_PROJECT_INST_DIR/$WM_PROJECT-$WM_FORK-$WM_PROJECT_VERSION +export WM_PROJECT_DIR="$projectDir" export WM_PROJECT_USER_DIR=$HOME/$WM_PROJECT/$USER-$WM_PROJECT_VERSION # Location of third-party software # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ #: ${WM_THIRD_PARTY_DIR=$WM_PROJECT_INST_DIR/ThirdParty-$WM_PROJECT_VERSION}; export WM_THIRD_PARTY_DIR +# Use external ThirdParty directory, e.g., from the standard installation of foam-extend. +#export WM_THIRD_PARTY_DIR=$HOME/$WM_PROJECT/$WM_PROJECT-$WM_FORK-4.1/ThirdParty export WM_THIRD_PARTY_DIR=$WM_PROJECT_DIR/ThirdParty # Enabling the usage of third-party software @@ -101,15 +101,15 @@ export WM_THIRD_PARTY_DIR=$WM_PROJECT_DIR/ThirdParty # # For AllMake.stage3 -export METIS_SYSTEM=1 -#export WM_THIRD_PARTY_USE_METIS_510=1 +#export METIS_SYSTEM=1 +export WM_THIRD_PARTY_USE_METIS_510=1 export WM_THIRD_PARTY_USE_PARMGRIDGEN_10=1 #export WM_THIRD_PARTY_USE_LIBCCMIO_261=1 export WM_THIRD_PARTY_USE_MESQUITE_230=1 -export SCOTCH_SYSTEM=1 -#export WM_THIRD_PARTY_USE_SCOTCH_604=1 -export PARMETIS_SYSTEM=1 -#export WM_THIRD_PARTY_USE_PARMETIS_403=1 +#export SCOTCH_SYSTEM=1 +export WM_THIRD_PARTY_USE_SCOTCH_604=1 +#export PARMETIS_SYSTEM=1 +export WM_THIRD_PARTY_USE_PARMETIS_403=1 export WM_THIRD_PARTY_USE_PYFOAM_069=1 export WM_THIRD_PARTY_USE_HWLOC_201=1 export WM_THIRD_PARTY_USE_ADIOS2_291=1 @@ -187,7 +187,7 @@ export WM_COMPILER_LIB_ARCH= : ${WM_PRECISION_OPTION:=DP}; export WM_PRECISION_OPTION # WM_LABEL_SIZE = 32 | 64 -: ${WM_LABEL_SIZE:=32}; export WM_LABEL_SIZE +: ${WM_LABEL_SIZE:=64}; export WM_LABEL_SIZE # WM_COMPILE_OPTION = Opt | Debug | Prof : ${WM_COMPILE_OPTION:=Opt}; export WM_COMPILE_OPTION @@ -601,7 +601,7 @@ fi # cleanup environment: # ~~~~~~~~~~~~~~~~~~~~ -unset cleanEnv cleanProg foamInstall foamOldDirs +unset cleanEnv cleanProg projectDir foamInstall foamOldDirs unset _foamSource # ----------------------------------------------------------------------------- diff --git a/src/decompositionMethods/decomposeReconstruct/decomposeTools/finiteVolume/decomposeMesh.C b/src/decompositionMethods/decomposeReconstruct/decomposeTools/finiteVolume/decomposeMesh.C index db73122ac..b803a9523 100644 --- a/src/decompositionMethods/decomposeReconstruct/decomposeTools/finiteVolume/decomposeMesh.C +++ b/src/decompositionMethods/decomposeReconstruct/decomposeTools/finiteVolume/decomposeMesh.C @@ -128,14 +128,6 @@ void Foam::domainDecomposition::decomposeMesh(const bool filterEmptyPatches) Pout<< "\nCalculating original mesh data" << endl; } - // Set references to the original mesh - const polyBoundaryMesh& patches = mesh_.boundaryMesh(); - - // Access all faces to grab the zones - const faceList& allFaces = mesh_.allFaces(); - const labelList& owner = mesh_.faceOwner(); - const labelList& neighbour = mesh_.faceNeighbour(); - // loop through the list of processor labels for the cell and add the // cell shape to the list of cells for the appropriate processor @@ -164,17 +156,32 @@ void Foam::domainDecomposition::decomposeMesh(const bool filterEmptyPatches) } // Convert linked lists into normal lists - forAll (procCellList, procI) + #pragma omp parallel for + for (label procI = 0; procI < procCellList.size(); ++procI) { procCellAddressing_[procI] = procCellList[procI]; } } + if (mesh_.time().writeFormat() == IOstream::COHERENT) + { + Info << "Skipping stones." << endl; + } + if (mesh_.time().writeFormat() != IOstream::COHERENT) + { if (debug) { Pout << "\nDistributing faces to processors" << endl; } + // Set references to the original mesh + const polyBoundaryMesh& patches = mesh_.boundaryMesh(); + + // Access all faces to grab the zones + const faceList& allFaces = mesh_.allFaces(); + const labelList& owner = mesh_.faceOwner(); + const labelList& neighbour = mesh_.faceNeighbour(); + // Loop through internal faces and decide which processor they belong to // First visit all internal faces. If cells at both sides belong to the // same processor, the face is an internal face. If they are different, @@ -1217,6 +1224,7 @@ void Foam::domainDecomposition::decomposeMesh(const bool filterEmptyPatches) globallySharedPoints_ = gSharedPoints.toc(); sort(globallySharedPoints_); } + } } diff --git a/src/decompositionMethods/decomposeReconstruct/decomposeTools/finiteVolume/domainDecomposition.C b/src/decompositionMethods/decomposeReconstruct/decomposeTools/finiteVolume/domainDecomposition.C index cdb9ea84a..0904f22bb 100644 --- a/src/decompositionMethods/decomposeReconstruct/decomposeTools/finiteVolume/domainDecomposition.C +++ b/src/decompositionMethods/decomposeReconstruct/decomposeTools/finiteVolume/domainDecomposition.C @@ -36,6 +36,9 @@ License #include "DynamicList.H" #include "globalMeshData.H" +#include "sliceMeshHelper.H" +#include "sliceWritePrimitives.H" + // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // defineTypeNameAndDebug(Foam::domainDecomposition, 0); @@ -87,6 +90,162 @@ Foam::domainDecomposition::~domainDecomposition() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // +Foam::autoPtr Foam::domainDecomposition::parallelMesh +( + const Time& processorDb, + const word& regionName +) +{ + // Create cell lookup + coherentCellAddressing_.resize(mesh_.nCells(), -1); + + label start = 0; + label gCellI = 0; + labelList partitionStarts( procCellAddressing_.size() + 1 ); + for (label procI = 0; procI < nProcs_; procI++) { + const labelList& curCellLabels = procCellAddressing_[procI]; + // Fill partitionStarts + partitionStarts[ procI ] = start; + start += curCellLabels.size(); + partitionStarts[ procI+1 ] = start; + // Fill the permutation from new to old cell numbering + forAll (curCellLabels, cellI) { + coherentCellAddressing_[curCellLabels[cellI]] = gCellI; + ++gCellI; + } + } + auto path = mesh_.pointsInstance()/mesh_.meshDir(); + sliceWritePrimitives + ( + "mesh", + path, + "partitionStarts", + partitionStarts.size(), + partitionStarts.cdata() + ); + + // Get complete owner-neighour addressing in the mesh + const labelList& own = mesh_.faceOwner(); + const labelList& nei = mesh_.faceNeighbour(); + + labelList procOwner(own.size(), -1); + labelList procNeighbour(nei.size(), -1); + coherentFaceAddressing_.resize(own.size(), 0); + faceList procFaces = mesh_.allFaces(); + + // Packaging mesh data to apply permutations + std::vector + < + std::tuple + > + owner_neighbour_face(nei.size()); + #pragma omp parallel for + for (Foam::label faceI = 0; faceI < nei.size(); ++faceI) + { + auto ownCellId = own[faceI]; + auto neiCellId = nei[faceI]; + auto owner = coherentCellAddressing_[ownCellId]; + auto neighbour = coherentCellAddressing_[neiCellId]; + auto face = procFaces[faceI]; + if ( owner > neighbour ) + { + std::swap(owner, neighbour); + face = face.reverseFace(); + } + owner_neighbour_face[faceI] = + std::make_tuple + ( + owner, + neighbour, + faceI, + face + ); + procNeighbour[faceI] = neighbour; + } + + labelList indices(procNeighbour.size()); + permutationOfSorted(indices.begin(), indices.end(), procNeighbour); + applyPermutation(owner_neighbour_face, indices); + + #pragma omp parallel for + for (Foam::label faceI = 0; faceI < owner_neighbour_face.size(); ++faceI) + { + procOwner[faceI] = std::get<0> + ( + std::move(owner_neighbour_face[faceI]) + ); + procNeighbour[faceI] = std::get<1> + ( + std::move(owner_neighbour_face[faceI]) + ); + coherentFaceAddressing_[faceI] = std::get<2> + ( + std::move(owner_neighbour_face[faceI]) + ); + procFaces[faceI] = std::get<3> + ( + std::move(owner_neighbour_face[faceI]) + ); + } + + // Filling owners for patch faces + #pragma omp parallel for + for (Foam::label faceI = nei.size(); faceI < own.size(); ++faceI) + { + auto cellId = own[faceI]; + procOwner[faceI] = coherentCellAddressing_[cellId]; + } + + // Filling face IDs for patch faces + std::iota + ( + coherentFaceAddressing_.begin() + procNeighbour.size(), + coherentFaceAddressing_.end(), + procNeighbour.size() + ); + + // Create processor mesh without a boundary + pointField procPoints = mesh_.allPoints(); + autoPtr procMeshPtr + ( + new fvMesh + ( + IOobject + ( + regionName, + mesh_.pointsInstance(), + processorDb + ), + xferMove(procPoints), + xferMove(procFaces), + xferMove(procOwner), + xferMove(procNeighbour), + false // Do not sync par + ) + ); + fvMesh& procMesh = procMeshPtr(); + + const polyPatchList& meshPatches = mesh_.boundaryMesh(); + List procPatches( meshPatches.size(), nullptr ); + forAll ( meshPatches, patchi ) + { + procPatches[patchi] = + meshPatches[patchi].clone + ( + procMesh.boundaryMesh() + ).ptr(); + } + procMesh.addFvPatches( procPatches, false ); + procMesh.write(); + + // TODO: Remove or insert as desired + //procMesh.checkMesh( true ); + // REMOVED: "Create processor boundary patches" + // Because will be identified through extended neighbour list while reading. + + return procMeshPtr; +} + Foam::autoPtr Foam::domainDecomposition::processorMesh ( const label procI, @@ -530,6 +689,21 @@ bool Foam::domainDecomposition::writeDecomposition() { Info<< "\nConstructing processor meshes" << endl; + if ( mesh_.time().writeFormat() == IOstream::COHERENT ) { + Time myDb + ( + Time::controlDictName, + mesh_.time().rootPath(), + mesh_.time().caseName(), + "system", + "constant", + true + ); + + parallelMesh( myDb, mesh_.polyMesh::name() ); + + } else { + // Make a lookup map for globally shared points Map