From 388ff4d91495d4f9b8c834350b0763de7131528c Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Mon, 23 Sep 2019 17:25:28 +0100 Subject: [PATCH] AMG support for overset meshes --- .../oversetAMGInterface/oversetAMGInterface.C | 305 ++++++++++++------ .../oversetAMGInterface/oversetAMGInterface.H | 14 +- .../oversetAMGInterfaceField.C | 37 ++- .../oversetFvPatch/oversetFvPatch.C | 17 +- .../oversetLduInterface/oversetLduInterface.H | 35 ++ .../oversetMesh/oversetMeshAddressing.C | 4 +- 6 files changed, 287 insertions(+), 125 deletions(-) diff --git a/src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.C b/src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.C index d991ff14b..d94754e9b 100644 --- a/src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.C +++ b/src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.C @@ -25,6 +25,8 @@ License #include "oversetAMGInterface.H" #include "oversetMesh.H" +#include "IPstream.H" +#include "OPstream.H" #include "addToRunTimeSelectionTable.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -60,64 +62,110 @@ void Foam::oversetAMGInterface::initMap() const << abort(FatalError); } - // Create send map: what I am sending to which processor + // Algorithm follows oversetMeshAddressing.C calcInterpolationMap + // with changes. Here, the information only exists on the acceptor + // side. HJ, 21/Sep/2019 + // Looping is performed over local acceptors - // There is no need to use a map: guess and resize - labelListList sendMap(Pstream::nProcs()); + const label nAcceptors = donorCellForAcceptor_.size(); - // Memory management - { - forAll (sendMap, procI) - { - sendMap[procI].setSize(donorCells_.size()); - } + // Record number of donors from any processor. This list will + // be reduced + labelListList nDonorsFromProcessorMap(Pstream::nProcs()); - // Count number of entries - labelList nSendMap(Pstream::nProcs(), 0); - - forAll (donorCellsProc_, dcpI) - { - // First index = proc - const label toProc = donorCellsProc_[dcpI]; - // Second index = fill size - // Content = donor cell index - sendMap[toProc][nSendMap[toProc]] = dcpI; - nSendMap[toProc]++; - } - - // Resize - forAll (sendMap, procI) - { - sendMap[procI].setSize(nSendMap[procI]); - } - } + // Collect donor cell index on the given processor + labelListList donorFromProcessor(Pstream::nProcs()); // Create construct map: what I am receiving from each processor labelListList constructMap(Pstream::nProcs()); - // Memory management + // Create list containing number of donors A processor is sending to other + // processors. + forAll (nDonorsFromProcessorMap, procI) { - forAll (constructMap, procI) + nDonorsFromProcessorMap[procI].setSize(Pstream::nProcs(), 0); + + donorFromProcessor[procI].setSize(nAcceptors); + constructMap[procI].setSize(nAcceptors); + } + + // Get the nDonorsFromProcessorMap for my processor + labelList& numberOfLocalDonorsFromProcs = + nDonorsFromProcessorMap[Pstream::myProcNo()]; + + labelList nConstructMap(Pstream::nProcs(), 0); + + forAll (donorCellForAcceptor_, accI) + { + // Record donor cell from donor proc + donorFromProcessor[donorProcForAcceptor_[accI]] + [numberOfLocalDonorsFromProcs[donorProcForAcceptor_[accI]]] = + donorProcIndexForAcceptor_[accI]; + + // Record that the processor is giving data + constructMap[donorProcForAcceptor_[accI]] + [numberOfLocalDonorsFromProcs[donorProcForAcceptor_[accI]]] = + accI; + + // Increment number of donors from this proc + numberOfLocalDonorsFromProcs[donorProcForAcceptor_[accI]]++; + } + + // Resize lists + forAll (donorFromProcessor, procI) + { + donorFromProcessor[procI].setSize + ( + numberOfLocalDonorsFromProcs[procI] + ); + + constructMap[procI].setSize + ( + numberOfLocalDonorsFromProcs[procI] + ); + } + + // Reduce the number of cells. This will allow the receiving + // processor to know it will receive the sendMap + Pstream::gatherList(nDonorsFromProcessorMap); + Pstream::scatterList(nDonorsFromProcessorMap); + + // Send the donor cell indices to other processor + // Lookup MY ROW + forAll (numberOfLocalDonorsFromProcs, nbrProcI) + { + // If there are entries, send them across + if (numberOfLocalDonorsFromProcs[nbrProcI] > 0) { - constructMap[procI].setSize(acceptorCells().size()); + OPstream toNbrProc + ( + Pstream::blocking, + nbrProcI + ); + + toNbrProc << donorFromProcessor[nbrProcI]; } + } - labelList nConstructMap(Pstream::nProcs(), 0); + // Create sendMap. All data will be received from nbr + // processors + labelListList sendMap(Pstream::nProcs()); - forAll (donorProcForAcceptor_, dpaI) + // Receive data from the other side. Loop through second + // index + forAll (nDonorsFromProcessorMap, nbrProcI) + { + if (nDonorsFromProcessorMap[nbrProcI][Pstream::myProcNo()] > 0) { - // First index = proc - const label fromProc = donorProcForAcceptor_[dpaI]; - // Second index = fill size - // Content = acceptor cell index - constructMap[fromProc][nConstructMap[fromProc]] = dpaI; - nConstructMap[fromProc]++; - } + // There is data to receive + IPstream fromNbrProc + ( + Pstream::blocking, + nbrProcI + ); - // Resize - forAll (constructMap, procI) - { - constructMap[procI].setSize(nConstructMap[procI]); + // Insert into sendMap for nbrProc + sendMap[nbrProcI] = labelList(fromNbrProc); } } @@ -127,6 +175,35 @@ void Foam::oversetAMGInterface::initMap() const sendMap, constructMap ); + + // // Collect donor cells. Use sendMap and unpack it + label nDonorCells = 0; + + forAll (sendMap, procI) + { + nDonorCells += sendMap[procI].size(); + } + + // Set donor cells and donor proc + donorCells_.setSize(nDonorCells); + donorCellsProc_.setSize(nDonorCells); + + // Reset counter + nDonorCells = 0; + + // Collect donor cells + forAll (sendMap, procI) + { + const labelList& curSend = sendMap[procI]; + + forAll (curSend, csI) + { + donorCells_[nDonorCells] = curSend[csI]; + donorCellsProc_[nDonorCells] = procI; + + nDonorCells++; + } + } } @@ -174,6 +251,17 @@ Foam::oversetAMGInterface::oversetAMGInterface // Coarsen local donor cells and renumber. Record local donor index // and communicate across + // Note + // For serial runs donors are calculated here + // For parallel runs, they will be over-ridden in initMap + // There are 2 ways of calculating donors + // - go to fine donors and pick up restrict addressing + // - take sendMap and unpack. The coarse donor will have the same value + // but the donor that provides data to multiple neighbouring processors + // will be recorded multiple times. However, its restrict index is + // the same. + // Reconsider. HJ, 23/Sep/2019 + // Get fine donors const labelList& fineDonors = fineOversetInterface_.donorCells(); @@ -182,32 +270,51 @@ Foam::oversetAMGInterface::oversetAMGInterface // Memory management { + // Note: donor collection must be done per-processor: a single + // donor may need to be passed to multiple acceptors on + // multiple processors. HJ, 23/Sep/2019 + // Into the hash table, using the coarse donor, record coarse processor // from fine donor. They must be the same - HashTable > coarseDonorProcMap; + List coarseDonorProcMap(Pstream::nProcs()); // For all fine donors, find and collect restrict addressing. // This will be the coarse donor forAll (fineDonors, fdI) { - coarseDonorProcMap.insert + coarseDonorProcMap[fineDonorsProc[fdI]].insert ( - localRestrictAddressing[fineDonors[fdI]], - fineDonorsProc[fdI] + localRestrictAddressing[fineDonors[fdI]] ); } // Collect coarse donor cells - donorCells_ = coarseDonorProcMap.sortedToc(); - // Pick up coarse donorCellsProc from the fine - donorCellsProc_.setSize(donorCells_.size()); + // Count donor cells + label nDonors = 0; - // This could be done with an iterator as well - // HJ, 13/Sep/2019 - forAll (donorCells_, dcI) + forAll (coarseDonorProcMap, procI) { - donorCellsProc_[dcI] = coarseDonorProcMap.find(donorCells_[dcI])(); + nDonors+= coarseDonorProcMap[procI].size(); + } + + donorCells_.setSize(nDonors); + donorCellsProc_.setSize(nDonors); + + // Reset counter + nDonors = 0; + + forAll (coarseDonorProcMap, procI) + { + const labelList& curDonors = coarseDonorProcMap[procI].sortedToc(); + + forAll (curDonors, cdI) + { + donorCells_[nDonors] = curDonors[cdI]; + donorCellsProc_[nDonors] = procI; + + nDonors++; + } } } @@ -243,6 +350,11 @@ Foam::oversetAMGInterface::oversetAMGInterface } // Communicate coarse donor + fineOversetInterface_.initInternalFieldTransfer + ( + Pstream::blocking, + acceptorCellDonorIndex + ); acceptorCellDonorIndex = fineOversetInterface_.internalFieldTransfer @@ -260,6 +372,12 @@ Foam::oversetAMGInterface::oversetAMGInterface ); // Communicate donor proc + fineOversetInterface_.initInternalFieldTransfer + ( + Pstream::blocking, + acceptorCellDonorProc + ); + acceptorCellDonorProc = fineOversetInterface_.internalFieldTransfer ( @@ -376,7 +494,7 @@ Foam::oversetAMGInterface::oversetAMGInterface curSlaveProcDonorId = acceptorCellDonorIndex[fringeFaceCells[ffI]]; - // Info<< "A: " << curMaster << tab << curMasterProc << tab << curSlave << tab << curSlaveProc << tab << curSlaveProcDonorId << endl; + // Pout<< "A: " << curMaster << tab << curMasterProc << tab << curSlave << tab << curSlaveProc << tab << curSlaveProcDonorId << endl; // Code in current master and slave curMaster += procOffset*curMasterProc; curSlave += procOffset*curSlaveProc; @@ -544,7 +662,7 @@ Foam::oversetAMGInterface::oversetAMGInterface curSlaveProc = acceptorCellDonorProc[ffI]; curSlaveProcDonorId = acceptorCellDonorIndex[ffI]; - // Info<< "B: " << curMaster << tab << curMasterProc << tab << curSlave << tab << curSlaveProc << tab << curSlaveProcDonorId << endl; + // Pout<< "B: " << curMaster << tab << curMasterProc << tab << curSlave << tab << curSlaveProc << tab << curSlaveProcDonorId << endl; // Code in current master and slave curMaster += procOffset*curMasterProc; curSlave += procOffset*curSlaveProc; @@ -698,7 +816,8 @@ Foam::oversetAMGInterface::oversetAMGInterface forAll (curNbrs, curNbrI) { // Record new coarse acceptor - faceCells_[nCoarseFaces] = curMaster; + faceCells_[nCoarseFaces] = + curMaster - procOffset*Pstream::myProcNo(); // Record donor for coarse acceptor (can be on different processor) donorCellForAcceptor_[nCoarseFaces] = curNbrs[curNbrI]; @@ -726,13 +845,18 @@ Foam::oversetAMGInterface::oversetAMGInterface } // Info<< "Finished level. Sizes " << nl - // << "faceCells_: " << faceCells_.size() + // << "faceCells_: " << faceCells_ // << " fineAddressing_: " << fineAddressing_.size() // << " restrictAddressing_: " << restrictAddressing_.size() // << " donorCells_: " << donorCells_.size() // << " donorCellForAcceptor_: " << donorCellForAcceptor_.size() // << " donorProcIndexForAcceptor_: " << donorProcIndexForAcceptor_.size() // << endl; + + if (Pstream::parRun()) + { + initMap(); + } } @@ -853,62 +977,63 @@ void Foam::oversetAMGInterface::initInternalFieldTransfer const Pstream::commsTypes commsType, const unallocLabelList& iF ) const -{} - - -Foam::tmp Foam::oversetAMGInterface::internalFieldTransfer -( - const Pstream::commsTypes commsType, - const unallocLabelList& iF -) const { // Repackage donor data to acceptors // Note: result is copied as internal field, re-scaled and passed across if (Pstream::parRun()) { - tmp tresult(new labelField(iF)); - labelList& result = tresult(); + labelTransferBuffer_ = iF; - map().distribute(result); - - return tresult; + map().distribute(labelTransferBuffer_); } else { - return tmp(new labelField(iF, donorCellForAcceptor_)); + labelTransferBuffer_ = labelField(iF, donorCellForAcceptor_); } } +Foam::tmp Foam::oversetAMGInterface::internalFieldTransfer +( + const Pstream::commsTypes commsType, + const unallocLabelList& +) const +{ + return tmp(new labelField(labelTransferBuffer_)); +} + + void Foam::oversetAMGInterface::initInternalFieldTransfer ( const Pstream::commsTypes commsType, const scalarField& iF ) const -{} +{ + // Repackage donor data to acceptors + // Note: result is copied as internal field, re-scaled and passed across + + // Repackage donor data to acceptors + // Note: result is copied as internal field, re-scaled and passed across + if (Pstream::parRun()) + { + fieldTransferBuffer_ = iF; + + map().distribute(fieldTransferBuffer_); + } + else + { + fieldTransferBuffer_ = scalarField(iF, donorCellForAcceptor_); + } +} Foam::tmp Foam::oversetAMGInterface::internalFieldTransfer ( const Pstream::commsTypes, - const scalarField& iF + const scalarField& ) const { - // Repackage donor data to acceptors - // Note: result is copied as internal field, re-scaled and passed across - if (Pstream::parRun()) - { - tmp tresult(new scalarField(iF)); - scalarList& result = tresult(); - - map().distribute(result); - - return tresult; - } - else - { - return tmp(new scalarField(iF, donorCellForAcceptor_)); - } + return tmp(new scalarField(fieldTransferBuffer_)); } diff --git a/src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.H b/src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.H index 864e88906..5784a5df9 100644 --- a/src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.H +++ b/src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.H @@ -64,7 +64,7 @@ class oversetAMGInterface //- Coarse interface size label interfaceSize_; - + // Local acceptor data // Local acceptor index is stored in faceCells_ @@ -79,18 +79,18 @@ class oversetAMGInterface //- Donor processor donor index for local acceptor labelList donorProcIndexForAcceptor_; - //- Local donor cells - labelList donorCells_; - - //- Donor cell processor (to which donor is sent) - labelList donorCellsProc_; - // Parallel communication //- Map-distribute comms tool mutable mapDistribute* mapPtr_; + //- Local donor cells + mutable labelList donorCells_; + + //- Donor cell processor (to which donor is sent) + mutable labelList donorCellsProc_; + // Private Member Functions diff --git a/src/overset/oversetMesh/oversetAMGInterfaceField/oversetAMGInterfaceField.C b/src/overset/oversetMesh/oversetAMGInterfaceField/oversetAMGInterfaceField.C index 0797ceddf..58d784160 100644 --- a/src/overset/oversetMesh/oversetAMGInterfaceField/oversetAMGInterfaceField.C +++ b/src/overset/oversetMesh/oversetAMGInterfaceField/oversetAMGInterfaceField.C @@ -58,12 +58,32 @@ Foam::oversetAMGInterfaceField::oversetAMGInterfaceField // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // void Foam::oversetAMGInterfaceField::initInterfaceMatrixUpdate +( + const scalarField& psiInternal, + scalarField&, + const lduMatrix& m, + const scalarField& coeffs, + const direction, + const Pstream::commsTypes commsType, + const bool +) const +{ + // Send psi to donors + oversetInterface_.initInternalFieldTransfer + ( + commsType, + psiInternal + ); +} + + +void Foam::oversetAMGInterfaceField::updateInterfaceMatrix ( const scalarField& psiInternal, scalarField& result, const lduMatrix& m, const scalarField& coeffs, - const direction cmpt, + const direction, const Pstream::commsTypes commsType, const bool switchToLhs ) const @@ -95,19 +115,4 @@ void Foam::oversetAMGInterfaceField::initInterfaceMatrixUpdate } -void Foam::oversetAMGInterfaceField::updateInterfaceMatrix -( - const scalarField&, - scalarField&, - const lduMatrix&, - const scalarField&, - const direction, - const Pstream::commsTypes commsType, - const bool switchToLhs -) const -{ - // Nothing to do -} - - // ************************************************************************* // diff --git a/src/overset/oversetMesh/oversetFvPatch/oversetFvPatch.C b/src/overset/oversetMesh/oversetFvPatch/oversetFvPatch.C index 4afd9dcb8..d342353da 100644 --- a/src/overset/oversetMesh/oversetFvPatch/oversetFvPatch.C +++ b/src/overset/oversetMesh/oversetFvPatch/oversetFvPatch.C @@ -159,7 +159,12 @@ void Foam::oversetFvPatch::initInternalFieldTransfer const Pstream::commsTypes commsType, const unallocLabelList& iF ) const -{} +{ + // Repackage donor data to acceptors + + // Note: result is copied as internal field, re-scaled and passed across + labelTransferBuffer_ = overset().acceptorMasterData(iF); +} Foam::tmp Foam::oversetFvPatch::internalFieldTransfer @@ -168,15 +173,7 @@ Foam::tmp Foam::oversetFvPatch::internalFieldTransfer const unallocLabelList& iF ) const { - // Repackage donor data to acceptors - - // Note: result is copied as internal field, re-scaled and passed across - tmp tresult(new labelField(iF)); - labelList& result = tresult(); - - result = overset().acceptorMasterData(result); - - return tresult; + return tmp(new labelField(labelTransferBuffer_)); } diff --git a/src/overset/oversetMesh/oversetLduInterface/oversetLduInterface.H b/src/overset/oversetMesh/oversetLduInterface/oversetLduInterface.H index d32f1239a..ccd72987e 100644 --- a/src/overset/oversetMesh/oversetLduInterface/oversetLduInterface.H +++ b/src/overset/oversetMesh/oversetLduInterface/oversetLduInterface.H @@ -40,6 +40,7 @@ SourceFiles #include "lduInterface.H" #include "primitiveFieldsFwd.H" +#include "crMatrix.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -57,6 +58,19 @@ class mapDistribute; class oversetLduInterface { +protected: + + // Protected data + + //- Transfer buffer + mutable labelField labelTransferBuffer_; + + //- Field transfer buffer + mutable scalarField fieldTransferBuffer_; + + //- crMatrix transfer buffer + mutable crMatrix crMatrixTransferBuffer_; + public: @@ -123,6 +137,27 @@ public: const Pstream::commsTypes commsType, const unallocLabelList& iF ) const = 0; + + + // Transfer buffer access + + //- Return contents of the label transfer buffer + const labelField& labelTransferBuffer() const + { + return labelTransferBuffer_; + } + + //- Return contents of the field transfer buffer + const scalarField& fieldTransferBuffer() const + { + return fieldTransferBuffer_; + } + + //- Return contents of the crMatrix transfer buffer + const crMatrix& crMatrixTransferBuffer() const + { + return crMatrixTransferBuffer_; + } }; diff --git a/src/overset/oversetMesh/oversetMesh/oversetMeshAddressing.C b/src/overset/oversetMesh/oversetMesh/oversetMeshAddressing.C index 05817faf6..19d6cecba 100644 --- a/src/overset/oversetMesh/oversetMesh/oversetMeshAddressing.C +++ b/src/overset/oversetMesh/oversetMesh/oversetMeshAddressing.C @@ -855,7 +855,7 @@ void Foam::oversetMesh::calcInterpolationMap() const // donor values I need to send to which processor: // - Loop through all regions and then through all donors for that region // - Mark all donors associated with donor/acceptor pair and see to which - // processor I need to send its data. Handle the possiblity that a + // processor I need to send its data. Handle the possibility that a // single donor may be used by multiple acceptors on the same processor // - While looping through donors, count how many donors I'm sending to // each processor @@ -957,7 +957,7 @@ void Foam::oversetMesh::calcInterpolationMap() const labelListList sendDataMap(Pstream::nProcs()); forAll (sendDataMap, procI) { - sendDataMap[procI] = sendMap[procI].toc(); + sendDataMap[procI] = sendMap[procI].sortedToc(); } // Gather/scatter number of donors going to each processor from each