AMG support for overset meshes

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
+    const label nAcceptors = donorCellForAcceptor_.size();
-    labelListList sendMap(Pstream::nProcs());
 
-    // Memory management
+    // Record number of donors from any processor.  This list will
-    {
+    // be reduced
-        forAll (sendMap, procI)
+    labelListList nDonorsFromProcessorMap(Pstream::nProcs());
-        {
-            sendMap[procI].setSize(donorCells_.size());
-        }
 
-        // Count number of entries
+    // Collect donor cell index on the given processor
-        labelList nSendMap(Pstream::nProcs(), 0);
+    labelListList donorFromProcessor(Pstream::nProcs());
-
-        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]);
-        }
-    }
 
     // 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);
-        {
+
-            constructMap[procI].setSize(acceptorCells().size());
+        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 (donorProcForAcceptor_, dpaI)
+    forAll (donorCellForAcceptor_, accI)
     {
-            // First index = proc
+        // Record donor cell from donor proc
-            const label fromProc = donorProcForAcceptor_[dpaI];
+        donorFromProcessor[donorProcForAcceptor_[accI]]
-            // Second index = fill size
+            [numberOfLocalDonorsFromProcs[donorProcForAcceptor_[accI]]] =
-            // Content = acceptor cell index
+            donorProcIndexForAcceptor_[accI];
-            constructMap[fromProc][nConstructMap[fromProc]] = dpaI;
+
-            nConstructMap[fromProc]++;
+        // 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
+    // Resize lists
-        forAll (constructMap, procI)
+    forAll (donorFromProcessor, procI)
     {
-            constructMap[procI].setSize(nConstructMap[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)
+        {
+            OPstream toNbrProc
+            (
+                Pstream::blocking,
+                nbrProcI
+            );
+
+            toNbrProc << donorFromProcessor[nbrProcI];
+        }
+    }
+
+    // Create sendMap.  All data will be received from nbr
+    // processors
+    labelListList sendMap(Pstream::nProcs());
+
+    // Receive data from the other side.  Loop through second
+    // index
+    forAll (nDonorsFromProcessorMap, nbrProcI)
+    {
+        if (nDonorsFromProcessorMap[nbrProcI][Pstream::myProcNo()] > 0)
+        {
+            // There is data to receive
+            IPstream fromNbrProc
+            (
+                Pstream::blocking,
+                nbrProcI
+            );
+
+            // 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<label, label, Hash<label> > coarseDonorProcMap;
+        List<labelHashSet> 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]],
+                localRestrictAddressing[fineDonors[fdI]]
-                fineDonorsProc[fdI]
             );
         }
 
         // Collect coarse donor cells
-        donorCells_ = coarseDonorProcMap.sortedToc();
 
-        // Pick up coarse donorCellsProc from the fine
+        // Count donor cells
-        donorCellsProc_.setSize(donorCells_.size());
+        label nDonors = 0;
 
-        // This could be done with an iterator as well
+        forAll (coarseDonorProcMap, procI)
-        // HJ, 13/Sep/2019
-        forAll (donorCells_, dcI)
         {
-            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::labelField> 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<labelField> tresult(new labelField(iF));
+        labelTransferBuffer_ = iF;
-        labelList& result = tresult();
 
-        map().distribute(result);
+        map().distribute(labelTransferBuffer_);
-
-        return tresult;
     }
     else
     {
-        return tmp<labelField>(new labelField(iF, donorCellForAcceptor_));
+        labelTransferBuffer_ = labelField(iF, donorCellForAcceptor_);
     }
 }
 
 
+Foam::tmp<Foam::labelField> Foam::oversetAMGInterface::internalFieldTransfer
+(
+    const Pstream::commsTypes commsType,
+    const unallocLabelList&
+) const
+{
+    return tmp<labelField>(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::scalarField> Foam::oversetAMGInterface::internalFieldTransfer
 (
     const Pstream::commsTypes,
-    const scalarField& iF
+    const scalarField&
 ) const
 {
-    // Repackage donor data to acceptors
+    return tmp<scalarField>(new scalarField(fieldTransferBuffer_));
-    // Note: result is copied as internal field, re-scaled and passed across
-    if (Pstream::parRun())
-    {
-        tmp<scalarField> tresult(new scalarField(iF));
-        scalarList& result = tresult();
-
-        map().distribute(result);
-
-        return tresult;
-    }
-    else
-    {
-        return tmp<scalarField>(new scalarField(iF, donorCellForAcceptor_));
-    }
 }
 
 

src/overset/oversetMesh/oversetAMGInterface/oversetAMGInterface.H

@@ -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
 

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
-}
-
-
 // ************************************************************************* //

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::labelField> Foam::oversetFvPatch::internalFieldTransfer
@@ -168,15 +173,7 @@ Foam::tmp<Foam::labelField> Foam::oversetFvPatch::internalFieldTransfer
     const unallocLabelList& iF
 ) const
 {
-    // Repackage donor data to acceptors
+    return tmp<labelField>(new labelField(labelTransferBuffer_));
-
-    // Note: result is copied as internal field, re-scaled and passed across
-    tmp<labelField> tresult(new labelField(iF));
-    labelList& result = tresult();
-    
-    result = overset().acceptorMasterData(result);
-
-    return tresult;
 }
 
 

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_;
+            }
 };
 
 

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