From 929a92bf7d1b8f90de0ad3e8e59e8e1c1fab74b6 Mon Sep 17 00:00:00 2001 From: mattijs Date: Wed, 8 May 2013 12:51:39 +0100 Subject: [PATCH] BUG: lduPrimitiveMesh: combining external patches in synchronised order --- .../GAMGAgglomerationTemplates.C | 14 +- .../GAMGProcAgglomeration.C | 150 +++++++++++++++++- .../GAMGProcAgglomeration.H | 6 +- .../processorGAMGInterface.H | 1 - .../meshes/lduMesh/lduPrimitiveMesh.C | 46 ++++-- .../polyMesh/globalMeshData/globalIndex.H | 52 +++++- .../globalMeshData/globalIndexTemplates.C | 40 ++--- 7 files changed, 253 insertions(+), 56 deletions(-) diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerationTemplates.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerationTemplates.C index e30a1e69cd..ce4fa0073a 100644 --- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerationTemplates.C +++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerationTemplates.C @@ -126,8 +126,7 @@ void Foam::GAMGAgglomeration::restrictField const List& procIDs = agglomProcIDs(coarseLevelIndex); const labelList& offsets = cellOffsets(coarseLevelIndex); - const globalIndex cellOffsetter(offsets); - cellOffsetter.gather(fineComm, procIDs, cf); + globalIndex::gather(offsets, fineComm, procIDs, cf); } } @@ -192,19 +191,10 @@ void Foam::GAMGAgglomeration::prolongField const List& procIDs = agglomProcIDs(coarseLevelIndex); const labelList& offsets = cellOffsets(coarseLevelIndex); - const globalIndex cellOffsetter(offsets); - label localSize = nCells_[levelIndex]; - Field allCf(localSize); - cellOffsetter.scatter - ( - coarseComm, - procIDs, - cf, - allCf - ); + globalIndex::scatter(offsets, coarseComm, procIDs, cf, allCf); forAll(fineToCoarse, i) { diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/GAMGProcAgglomeration/GAMGProcAgglomeration.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/GAMGProcAgglomeration/GAMGProcAgglomeration.C index f825ca769d..d3a83dac79 100644 --- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/GAMGProcAgglomeration/GAMGProcAgglomeration.C +++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/GAMGProcAgglomeration/GAMGProcAgglomeration.C @@ -25,6 +25,7 @@ License #include "GAMGProcAgglomeration.H" #include "GAMGAgglomeration.H" +#include "lduMesh.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -114,6 +115,152 @@ void Foam::GAMGProcAgglomeration::printStats } +Foam::labelListList Foam::GAMGProcAgglomeration::globalCellCells +( + const lduMesh& mesh +) +{ + const lduAddressing& addr = mesh.lduAddr(); + lduInterfacePtrsList interfaces = mesh.interfaces(); + + const label myProcID = Pstream::myProcNo(mesh.comm()); + + globalIndex globalNumbering + ( + addr.size(), + Pstream::msgType(), + mesh.comm(), + Pstream::parRun() + ); + + labelList globalIndices(addr.size()); + forAll(globalIndices, cellI) + { + globalIndices[cellI] = globalNumbering.toGlobal(myProcID, cellI); + } + + + // Get the interface cells + PtrList nbrGlobalCells(interfaces.size()); + { + // Initialise transfer of restrict addressing on the interface + forAll(interfaces, inti) + { + if (interfaces.set(inti)) + { + interfaces[inti].initInternalFieldTransfer + ( + Pstream::nonBlocking, + globalIndices + ); + } + } + + if (Pstream::parRun()) + { + Pstream::waitRequests(); + } + + forAll(interfaces, inti) + { + if (interfaces.set(inti)) + { + nbrGlobalCells.set + ( + inti, + new labelList + ( + interfaces[inti].internalFieldTransfer + ( + Pstream::nonBlocking, + globalIndices + ) + ) + ); + } + } + } + + + // Scan the neighbour list to find out how many times the cell + // appears as a neighbour of the face. Done this way to avoid guessing + // and resizing list + labelList nNbrs(addr.size(), 1); + + const labelUList& nbr = addr.upperAddr(); + const labelUList& own = addr.lowerAddr(); + + { + forAll(nbr, faceI) + { + nNbrs[nbr[faceI]]++; + nNbrs[own[faceI]]++; + } + + forAll(interfaces, inti) + { + if (interfaces.set(inti)) + { + const labelUList& faceCells = interfaces[inti].faceCells(); + + forAll(faceCells, i) + { + nNbrs[faceCells[i]]++; + } + } + } + } + + + // Create cell-cells addressing + labelListList cellCells(addr.size()); + + forAll(cellCells, cellI) + { + cellCells[cellI].setSize(nNbrs[cellI], -1); + } + + // Reset the list of number of neighbours to zero + nNbrs = 0; + + // Scatter the neighbour faces + forAll(nbr, faceI) + { + label c0 = own[faceI]; + label c1 = nbr[faceI]; + + cellCells[c0][nNbrs[c0]++] = globalIndices[c1]; + cellCells[c1][nNbrs[c1]++] = globalIndices[c0]; + } + forAll(interfaces, inti) + { + if (interfaces.set(inti)) + { + const labelUList& faceCells = interfaces[inti].faceCells(); + + forAll(faceCells, i) + { + label c0 = faceCells[i]; + cellCells[c0][nNbrs[c0]++] = nbrGlobalCells[inti][i]; + } + } + } + + forAll(cellCells, cellI) + { + Foam::stableSort(cellCells[cellI]); + } + + // Replace the initial element (always -1) with the local cell + forAll(cellCells, cellI) + { + cellCells[cellI][0] = globalIndices[cellI]; + } + + return cellCells; +} + + bool Foam::GAMGProcAgglomeration::agglomerate ( const label fineLevelIndex, @@ -243,7 +390,4 @@ Foam::GAMGProcAgglomeration::~GAMGProcAgglomeration() {} -// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // - - // ************************************************************************* // diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/GAMGProcAgglomeration/GAMGProcAgglomeration.H b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/GAMGProcAgglomeration/GAMGProcAgglomeration.H index cf27f44bc0..c3acc01521 100644 --- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/GAMGProcAgglomeration/GAMGProcAgglomeration.H +++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/GAMGProcAgglomeration/GAMGProcAgglomeration.H @@ -44,6 +44,7 @@ namespace Foam { class GAMGAgglomeration; +class lduMesh; /*---------------------------------------------------------------------------*\ Class GAMGProcAgglomeration Declaration @@ -74,6 +75,9 @@ protected: const label procAgglomComm ); + //- Debug: calculate global cell-cells + static labelListList globalCellCells(const lduMesh&); + private: // Private data @@ -139,7 +143,7 @@ public: // Member Functions - //- Modify agglomeration. Return true if modified + //- Modify agglomeration. Return true if modified virtual bool agglomerate() = 0; }; diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/interfaces/processorGAMGInterface/processorGAMGInterface.H b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/interfaces/processorGAMGInterface/processorGAMGInterface.H index df868544d0..3d0c9d28c3 100644 --- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/interfaces/processorGAMGInterface/processorGAMGInterface.H +++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/interfaces/processorGAMGInterface/processorGAMGInterface.H @@ -29,7 +29,6 @@ Description SourceFiles processorGAMGInterface.C - processorGAMGInterfaceTemplates.C \*---------------------------------------------------------------------------*/ diff --git a/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.C b/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.C index 68069c88cf..070edd4b69 100644 --- a/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.C +++ b/src/OpenFOAM/meshes/lduMesh/lduPrimitiveMesh.C @@ -360,6 +360,7 @@ Foam::lduPrimitiveMesh::lduPrimitiveMesh + procMesh.lduAddr().size(); } + // Faces initially get added in order, sorted later labelList internalFaceOffsets(nMeshes+1); internalFaceOffsets[0] = 0; @@ -818,9 +819,37 @@ Foam::lduPrimitiveMesh::lduPrimitiveMesh { const labelPairList& elems = iter(); - // Sort processors in increasing order - labelList order(identity(elems.size())); - Foam::sort(order, procLess(elems)); + // Sort processors in increasing order so both sides walk through in + // same order. + labelPairList procPairs(elems.size()); + forAll(elems, i) + { + const labelPair& elem = elems[i]; + label procMeshI = elem[0]; + label interfaceI = elem[1]; + const lduInterfacePtrsList interfaces = mesh + ( + myMesh, + otherMeshes, + procMeshI + ).interfaces(); + + const processorLduInterface& pldui = + refCast + ( + interfaces[interfaceI] + ); + label myProcNo = pldui.myProcNo(); + label nbrProcNo = pldui.neighbProcNo(); + procPairs[i] = labelPair + ( + min(myProcNo, nbrProcNo), + max(myProcNo, nbrProcNo) + ); + } + + labelList order; + sortedOrder(procPairs, order); // Count label n = 0; @@ -837,11 +866,6 @@ Foam::lduPrimitiveMesh::lduPrimitiveMesh procMeshI ).interfaces(); - //Pout<< " adding interface:" << " procMeshI:" << procMeshI - // << " interface:" << interfaceI - // << " size:" << interfaces[interfaceI].faceCells().size() - // << endl; - n += interfaces[interfaceI].faceCells().size(); } @@ -868,6 +892,7 @@ Foam::lduPrimitiveMesh::lduPrimitiveMesh const labelUList& l = interfaces[interfaceI].faceCells(); bfMap.setSize(l.size()); + forAll(l, faceI) { allFaceCells[n] = cellOffsets[procMeshI]+l[faceI]; @@ -947,7 +972,10 @@ Foam::lduPrimitiveMesh::lduPrimitiveMesh patchSchedule_ = nonBlockingSchedule(interfaces_); - checkUpperTriangular(cellOffsets.last(), lowerAddr_, upperAddr_); + if (debug) + { + checkUpperTriangular(cellOffsets.last(), lowerAddr_, upperAddr_); + } } diff --git a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndex.H b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndex.H index 4e92401489..80263b91d3 100644 --- a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndex.H +++ b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndex.H @@ -150,6 +150,19 @@ public: // Other + //- Collect data in processor order on master (== procIDs[0]). + // Needs offsets only on master. + template + static void gather + ( + const labelUList& offsets, + const label comm, + const labelList& procIDs, + const UList& fld, + List& allFld, + const int tag = UPstream::msgType() + ); + //- Collect data in processor order on master (== procIDs[0]). // Needs offsets only on master. template @@ -160,7 +173,22 @@ public: const UList& fld, List& allFld, const int tag = UPstream::msgType() - ) const; + ) const + { + gather(offsets_, comm, procIDs, fld, allFld, tag); + } + + //- Inplace collect data in processor order on master + // (== procIDs[0]). Needs offsets only on master. + template + static void gather + ( + const labelUList& offsets, + const label comm, + const labelList& procIDs, + List& fld, + const int tag = UPstream::msgType() + ); //- Inplace collect data in processor order on master // (== procIDs[0]). Needs offsets only on master. @@ -171,7 +199,22 @@ public: const labelList& procIDs, List& fld, const int tag = UPstream::msgType() - ) const; + ) const + { + gather(offsets_, comm, procIDs, fld, tag); + } + + //- Distribute data in processor order. Requires fld to be sized! + template + static void scatter + ( + const labelUList& offsets, + const label comm, + const labelList& procIDs, + const UList& allFld, + UList& fld, + const int tag = UPstream::msgType() + ); //- Distribute data in processor order. Requires fld to be sized! template @@ -182,7 +225,10 @@ public: const UList& allFld, UList& fld, const int tag = UPstream::msgType() - ) const; + ) const + { + scatter(offsets_, comm, procIDs, allFld, fld, tag); + } // IOstream Operators diff --git a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndexTemplates.C b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndexTemplates.C index 2beb2cfe7a..73d2ea25a8 100644 --- a/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndexTemplates.C +++ b/src/OpenFOAM/meshes/polyMesh/globalMeshData/globalIndexTemplates.C @@ -30,28 +30,24 @@ License template void Foam::globalIndex::gather ( + const labelUList& offsets, const label comm, const labelList& procIDs, const UList& fld, List& allFld, const int tag -) const +) { if (Pstream::myProcNo(comm) == procIDs[0]) { - allFld.setSize(size()); + allFld.setSize(offsets.last()); // Assign my local data SubList(allFld, fld.size(), 0).assign(fld); for (label i = 1; i < procIDs.size(); i++) { - SubList procSlot - ( - allFld, - offsets_[i+1]-offsets_[i], - offsets_[i] - ); + SubList procSlot(allFld, offsets[i+1]-offsets[i], offsets[i]); if (contiguous()) { @@ -112,27 +108,23 @@ void Foam::globalIndex::gather template void Foam::globalIndex::gather ( + const labelUList& offsets, const label comm, const labelList& procIDs, List& fld, const int tag -) const +) { if (Pstream::myProcNo(comm) == procIDs[0]) { - List allFld(size()); + List allFld(offsets.last()); // Assign my local data SubList(allFld, fld.size(), 0).assign(fld); for (label i = 1; i < procIDs.size(); i++) { - SubList procSlot - ( - allFld, - offsets_[i+1]-offsets_[i], - offsets_[i] - ); + SubList procSlot(allFld, offsets[i+1]-offsets[i], offsets[i]); if (contiguous()) { @@ -195,31 +187,25 @@ void Foam::globalIndex::gather template void Foam::globalIndex::scatter ( + const labelUList& offsets, const label comm, const labelList& procIDs, const UList& allFld, UList& fld, const int tag -) const +) { if (Pstream::myProcNo(comm) == procIDs[0]) { - fld.assign - ( - SubList - ( - allFld, - offsets_[1]-offsets_[0] - ) - ); + fld.assign(SubList(allFld, offsets[1]-offsets[0])); for (label i = 1; i < procIDs.size(); i++) { const SubList procSlot ( allFld, - offsets_[i+1]-offsets_[i], - offsets_[i] + offsets[i+1]-offsets[i], + offsets[i] ); if (contiguous())