Merge remote-tracking branch 'origin/feature/procAgglom'

etc/controlDict

@@ -135,7 +135,7 @@ DebugSwitches
     FDIC                0;
     FaceCellWave        0;
     GAMG                0;
-    GAMGAgglomeration   0;
+    GAMGAgglomeration   1;
     GAMGInterface       0;
     GAMGInterfaceField  0;
     Gamma               0;

src/OpenFOAM/Make/files

@@ -331,10 +331,8 @@ eagerGAMGProcAgglomeration = $(GAMGProcAgglomerations)/eagerGAMGProcAgglomeratio
 $(eagerGAMGProcAgglomeration)/eagerGAMGProcAgglomeration.C
 noneGAMGProcAgglomeration = $(GAMGProcAgglomerations)/noneGAMGProcAgglomeration
 $(noneGAMGProcAgglomeration)/noneGAMGProcAgglomeration.C
-/*
-cellFaceRatioGAMGProcAgglomeration = $(GAMGProcAgglomerations)/cellFaceRatioGAMGProcAgglomeration
-$(cellFaceRatioGAMGProcAgglomeration)/cellFaceRatioGAMGProcAgglomeration.C
-*/
+procFacesGAMGProcAgglomeration = $(GAMGProcAgglomerations)/procFacesGAMGProcAgglomeration
+$(procFacesGAMGProcAgglomeration)/procFacesGAMGProcAgglomeration.C
 
 
 meshes/lduMesh/lduMesh.C

src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomerateLduAddressing.C

@@ -351,7 +351,7 @@ void Foam::GAMGAgglomeration::agglomerateLduAddressing
     );
 
 
-    if (debug)
+    if (debug & 2)
     {
         Pout<< "GAMGAgglomeration :"
             << " agglomerated level " << fineLevelIndex

src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/GAMGAgglomeration/GAMGAgglomeration.C

@@ -27,9 +27,9 @@ License
 #include "lduMesh.H"
 #include "lduMatrix.H"
 #include "Time.H"
-#include "dlLibraryTable.H"
 #include "GAMGInterface.H"
 #include "GAMGProcAgglomeration.H"
+#include "IOmanip.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -67,23 +67,94 @@ void Foam::GAMGAgglomeration::compactLevels(const label nCreatedLevels)
         procBoundaryFaceMap_.setSize(nCreatedLevels);
 
         procAgglomeratorPtr_().agglomerate();
-    }
 
-    if (debug)
-    {
-        for (label levelI = 0; levelI <= size(); levelI++)
+
+        if (debug)
         {
-            if (hasMeshLevel(levelI))
-            {
-                const lduMesh& fineMesh = meshLevel(levelI);
-                Pout<< "Level " << levelI << " fine mesh:"<< nl;
-                Pout<< fineMesh.info() << endl;
-            }
-            else
+
+            Info<< "GAMGAgglomeration:" << nl
+                << "    local agglomerator     : " << type() << nl
+                << "    processor agglomerator : "
+                << procAgglomeratorPtr_().type() << nl
+                << nl;
+
+            Info<< setw(40) << "nCells"
+                << setw(24) << "nInterfaces"
+                << setw(24) << "Ratio" << nl
+                << setw(8) << "Level"
+                << setw(8) << "nProcs"
+                << "        "
+                << setw(8) << "avg"
+                << setw(8) << "max"
+                << "        "
+                << setw(8) << "avg"
+                << setw(8) << "max"
+                << "            " << setw(4) << "avg"
+                << "    " << setw(4) << "max"
+                << nl
+                << setw(8) << "-----"
+                << setw(8) << "------"
+                << "        "
+                << setw(8) << "---"
+                << setw(8) << "---"
+                << "        "
+                << setw(8) << "---"
+                << setw(8) << "---"
+                << "            " << setw(4) << "---"
+                << "    " << setw(4) << "---"
+                << nl;
+
+            for (label levelI = 0; levelI <= size(); levelI++)
             {
-                Pout<< "Level " << levelI << " has no fine mesh:" << nl
-                    << endl;
+                label nProcs = 0;
+                label nCells = 0;
+                label nInterfaces = 0;
+                label nIntFaces = 0;
+                scalar ratio = 0.0;
+
+                if (hasMeshLevel(levelI))
+                {
+                    nProcs = 1;
+
+                    const lduMesh& fineMesh = meshLevel(levelI);
+                    nCells = fineMesh.lduAddr().size();
+
+                    const lduInterfacePtrsList interfaces =
+                        fineMesh.interfaces();
+                    forAll(interfaces, i)
+                    {
+                        if (interfaces.set(i))
+                        {
+                            nInterfaces++;
+                            nIntFaces += interfaces[i].faceCells().size();
+                        }
+                    }
+                    ratio = scalar(nIntFaces)/nCells;
+                }
+
+                label totNprocs = returnReduce(nProcs, sumOp<label>());
+
+                label maxNCells = returnReduce(nCells, maxOp<label>());
+                label totNCells = returnReduce(nCells, sumOp<label>());
+
+                label maxNInt = returnReduce(nInterfaces, maxOp<label>());
+                label totNInt = returnReduce(nInterfaces, sumOp<label>());
+
+                scalar maxRatio = returnReduce(ratio, maxOp<scalar>());
+                scalar totRatio = returnReduce(ratio, sumOp<scalar>());
+
+                Info<< setw(8) << levelI
+                    << setw(8) << totNprocs
+                    << setw(16) << totNCells/totNprocs
+                    << setw(8) << maxNCells
+                    << setw(16) << totNInt/totNprocs
+                    << setw(8) << maxNInt
+                    << "        "
+                    << setw(8) << setprecision(4) << totRatio/totNprocs
+                    << setw(8) << setprecision(4) << maxRatio
+                    << nl;
             }
+            Info<< endl;
         }
     }
 }

src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGAgglomerations/pairGAMGAgglomeration/pairGAMGAgglomeration.H

@@ -61,14 +61,6 @@ protected:
 
     // Protected Member Functions
 
-        //- Calculate and return agglomeration of given level
-        tmp<labelField> agglomerate
-        (
-            label& nCoarseCells,
-            const lduAddressing& fineMatrixAddressing,
-            const scalarField& faceWeights
-        );
-
         //- Agglomerate all levels starting from the given face weights
         void agglomerate
         (
@@ -97,6 +89,15 @@ public:
             const lduMesh& mesh,
             const dictionary& controlDict
         );
+
+        //- Calculate and return agglomeration
+        static tmp<labelField> agglomerate
+        (
+            label& nCoarseCells,
+            const lduAddressing& fineMatrixAddressing,
+            const scalarField& faceWeights
+        );
+
 };
 
 

src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/procFacesGAMGProcAgglomeration/procFacesGAMGProcAgglomeration.C

+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "procFacesGAMGProcAgglomeration.H"
+#include "addToRunTimeSelectionTable.H"
+#include "GAMGAgglomeration.H"
+#include "Random.H"
+#include "lduMesh.H"
+#include "processorLduInterface.H"
+#include "processorGAMGInterface.H"
+#include "pairGAMGAgglomeration.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(procFacesGAMGProcAgglomeration, 0);
+
+    addToRunTimeSelectionTable
+    (
+        GAMGProcAgglomeration,
+        procFacesGAMGProcAgglomeration,
+        GAMGAgglomeration
+    );
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+// Create single cell mesh
+Foam::autoPtr<Foam::lduPrimitiveMesh>
+Foam::procFacesGAMGProcAgglomeration::singleCellMesh
+(
+    const label singleCellMeshComm,
+    const lduMesh& mesh,
+    scalarField& faceWeights
+) const
+{
+    // Count number of faces per processor
+    List<Map<label> > procFaces(UPstream::nProcs(mesh.comm()));
+    Map<label>& myNeighbours = procFaces[UPstream::myProcNo(mesh.comm())];
+
+    {
+        const lduInterfacePtrsList interfaces(mesh.interfaces());
+        forAll(interfaces, intI)
+        {
+            if (interfaces.set(intI))
+            {
+                const processorLduInterface& pp =
+                    refCast<const processorLduInterface>
+                    (
+                        interfaces[intI]
+                    );
+                label size = interfaces[intI].faceCells().size();
+                myNeighbours.insert(pp.neighbProcNo(), size);
+            }
+        }
+    }
+
+    Pstream::gatherList(procFaces, Pstream::msgType(), mesh.comm());
+    Pstream::scatterList(procFaces, Pstream::msgType(), mesh.comm());
+
+    autoPtr<lduPrimitiveMesh> singleCellMeshPtr;
+
+    if (Pstream::master(mesh.comm()))
+    {
+        // I am master
+        label nCells = Pstream::nProcs(mesh.comm());
+
+        DynamicList<label> l(3*nCells);
+        DynamicList<label> u(3*nCells);
+        DynamicList<scalar> weight(3*nCells);
+
+        DynamicList<label> nbrs;
+        DynamicList<scalar> weights;
+
+        forAll(procFaces, procI)
+        {
+            const Map<label>& neighbours = procFaces[procI];
+
+            // Add all the higher processors
+            nbrs.clear();
+            weights.clear();
+            forAllConstIter(Map<label>, neighbours, iter)
+            {
+                if (iter.key() > procI)
+                {
+                    nbrs.append(iter.key());
+                    weights.append(iter());
+                }
+                sort(nbrs);
+                forAll(nbrs, i)
+                {
+                    l.append(procI);
+                    u.append(nbrs[i]);
+                    weight.append(weights[i]);
+                }
+            }
+        }
+
+        faceWeights.transfer(weight);
+
+        PtrList<const lduInterface> primitiveInterfaces(0);
+        const lduSchedule ps(0);
+
+        singleCellMeshPtr.reset
+        (
+            new lduPrimitiveMesh
+            (
+                nCells,
+                l,
+                u,
+                primitiveInterfaces,
+                ps,
+                singleCellMeshComm
+            )
+        );
+    }
+    return singleCellMeshPtr;
+}
+
+
+Foam::tmp<Foam::labelField>
+Foam::procFacesGAMGProcAgglomeration::processorAgglomeration
+(
+    const lduMesh& mesh
+) const
+{
+    label singleCellMeshComm = UPstream::allocateCommunicator
+    (
+        mesh.comm(),
+        labelList(1, 0)            // only processor 0
+    );
+
+    scalarField faceWeights;
+    autoPtr<lduPrimitiveMesh> singleCellMeshPtr
+    (
+        singleCellMesh
+        (
+            singleCellMeshComm,
+            mesh,
+            faceWeights
+        )
+    );
+
+    tmp<labelField> tfineToCoarse(new labelField(0));
+    labelField& fineToCoarse = tfineToCoarse();
+
+    if (singleCellMeshPtr.valid())
+    {
+        // On master call the agglomerator
+        const lduPrimitiveMesh& singleCellMesh = singleCellMeshPtr();
+
+        label nCoarseProcs;
+        fineToCoarse = pairGAMGAgglomeration::agglomerate
+        (
+            nCoarseProcs,
+            singleCellMesh,
+            faceWeights
+        );
+
+        labelList coarseToMaster(nCoarseProcs, labelMax);
+        forAll(fineToCoarse, cellI)
+        {
+            label coarseI = fineToCoarse[cellI];
+            coarseToMaster[coarseI] = min(coarseToMaster[coarseI], cellI);
+        }
+
+        // Sort according to master and redo restriction
+        labelList newToOld;
+        sortedOrder(coarseToMaster, newToOld);
+        labelList oldToNew(invert(newToOld.size(), newToOld));
+
+        fineToCoarse = UIndirectList<label>(oldToNew, fineToCoarse)();
+    }
+
+    Pstream::scatter(fineToCoarse, Pstream::msgType(), mesh.comm());
+    UPstream::freeCommunicator(singleCellMeshComm);
+
+    return tfineToCoarse;
+}
+
+
+bool Foam::procFacesGAMGProcAgglomeration::doProcessorAgglomeration
+(
+    const lduMesh& mesh
+) const
+{
+    // Check the need for further agglomeration on all processors
+    bool doAgg = mesh.lduAddr().size() < nAgglomeratingCells_;
+    mesh.reduce(doAgg, orOp<bool>());
+    return doAgg;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::procFacesGAMGProcAgglomeration::procFacesGAMGProcAgglomeration
+(
+    GAMGAgglomeration& agglom,
+    const dictionary& controlDict
+)
+:
+    GAMGProcAgglomeration(agglom, controlDict),
+    nAgglomeratingCells_(readLabel(controlDict.lookup("nAgglomeratingCells")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::procFacesGAMGProcAgglomeration::~procFacesGAMGProcAgglomeration()
+{
+    forAllReverse(comms_, i)
+    {
+        if (comms_[i] != -1)
+        {
+            UPstream::freeCommunicator(comms_[i]);
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+bool Foam::procFacesGAMGProcAgglomeration::agglomerate()
+{
+    if (debug)
+    {
+        Pout<< nl << "Starting mesh overview" << endl;
+        printStats(Pout, agglom_);
+    }
+
+    if (agglom_.size() >= 1)
+    {
+        Random rndGen(0);
+
+        for
+        (
+            label fineLevelIndex = 2;
+            fineLevelIndex < agglom_.size();
+            fineLevelIndex++
+        )
+        {
+            if (agglom_.hasMeshLevel(fineLevelIndex))
+            {
+                // Get the fine mesh
+                const lduMesh& levelMesh = agglom_.meshLevel(fineLevelIndex);
+
+                label levelComm = levelMesh.comm();
+                label nProcs = UPstream::nProcs(levelComm);
+
+                if (nProcs > 1 && doProcessorAgglomeration(levelMesh))
+                {
+                    tmp<labelField> tprocAgglomMap
+                    (
+                        processorAgglomeration(levelMesh)
+                    );
+                    const labelField& procAgglomMap = tprocAgglomMap();
+
+                    // Master processor
+                    labelList masterProcs;
+                    // Local processors that agglomerate. agglomProcIDs[0] is in
+                    // masterProc.
+                    List<int> agglomProcIDs;
+                    GAMGAgglomeration::calculateRegionMaster
+                    (
+                        levelComm,
+                        procAgglomMap,
+                        masterProcs,
+                        agglomProcIDs
+                    );
+
+                    // Allocate a communicator for the processor-agglomerated
+                    // matrix
+                    comms_.append
+                    (
+                        UPstream::allocateCommunicator
+                        (
+                            levelComm,
+                            masterProcs
+                        )
+                    );
+
+
+                    // Use procesor agglomeration maps to do the actual
+                    // collecting.
+                    GAMGProcAgglomeration::agglomerate
+                    (
+                        fineLevelIndex,
+                        procAgglomMap,
+                        masterProcs,
+                        agglomProcIDs,
+                        comms_.last()
+                    );
+                }
+            }
+        }
+    }
+
+    // Print a bit
+    if (debug)
+    {
+        Pout<< nl << "Agglomerated mesh overview" << endl;
+        printStats(Pout, agglom_);
+    }
+
+    return true;
+}
+
+
+// ************************************************************************* //

src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGProcAgglomerations/procFacesGAMGProcAgglomeration/procFacesGAMGProcAgglomeration.H

+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::procFacesGAMGProcAgglomeration
+
+Description
+    Processor agglomeration of GAMGAgglomerations. Needs nAgglomeratingCells
+    which is when to start agglomerating processors. Processors get agglomerated
+    by constructing a single cell mesh for each processor with each
+    processor interface a face. This then gets agglomerated using
+    the pairGAMGAgglomeration algorithm with the number of faces
+    on the original processor interface as face weight.