Merge branch 'feature-kahip-decompose' into 'develop'

add support for KaHIP decomposition

See merge request Development/OpenFOAM-plus!154

applications/test/checkDecomposePar/Make/files

@@ -0,0 +1,3 @@
+Test-checkDecomposePar.C
+
+EXE = $(FOAM_APPBIN)/Test-checkDecomposePar

applications/test/checkDecomposePar/Make/options

@@ -0,0 +1,10 @@
+EXE_INC = \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/regionModels/regionModel/lnInclude
+
+EXE_LIBS = \
+    -ldecompositionMethods \
+    -lregionModels

applications/test/checkDecomposePar/Test-checkDecomposePar.C

@@ -0,0 +1,193 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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/>.
+
+Application
+    checkDecomposePar
+
+Group
+    grpParallelUtilities
+
+Description
+    Check decomposition from kaffpa (KaHIP) output.
+    foamToMetisGraph was likely used for producing the kaffpa input.
+
+\*---------------------------------------------------------------------------*/
+
+#include "OSspecific.H"
+#include "fvCFD.H"
+#include "cpuTime.H"
+#include "IFstream.H"
+#include "regionProperties.H"
+#include "decompositionInformation.H"
+#include "decompositionModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+int main(int argc, char *argv[])
+{
+    argList::addNote
+    (
+        "Check decomposition from kaffpa (KaHIP) output"
+    );
+
+    argList::noParallel();
+    argList::noBanner();
+
+    #include "addRegionOption.H"
+    argList::addBoolOption
+    (
+        "allRegions",
+        "operate on all regions in regionProperties"
+    );
+    argList::addBoolOption
+    (
+        "verbose",
+        "more information about decomposition"
+    );
+
+    argList::validArgs.append("kaffpa-output-file");
+
+    // Include explicit constant options, have zero from time range
+    timeSelector::addOptions(true, false);
+
+    #include "setRootCase.H"
+
+    const fileName decompFile = args[1];
+
+    const bool region     = args.optionFound("region");
+    const bool allRegions = args.optionFound("allRegions");
+    const bool verbose    = args.optionFound("verbose");
+
+    // Set time from database
+    #include "createTime.H"
+    // Allow override of time
+    instantList times = timeSelector::selectIfPresent(runTime, args);
+
+    // Allow override of decomposeParDict location
+    fileName decompDictFile;
+    args.optionReadIfPresent("decomposeParDict", decompDictFile);
+
+    wordList regionNames;
+    wordList regionDirs;
+    if (allRegions)
+    {
+        Info<< "Decomposing all regions in regionProperties" << nl << endl;
+        regionProperties rp(runTime);
+        forAllConstIters(rp, iter)
+        {
+            const wordList& regions = iter();
+            forAll(regions, i)
+            {
+                if (findIndex(regionNames, regions[i]) == -1)
+                {
+                    regionNames.append(regions[i]);
+                }
+            }
+        }
+        regionDirs = regionNames;
+    }
+    else
+    {
+        word regionName;
+        if (args.optionReadIfPresent("region", regionName))
+        {
+            regionNames = wordList(1, regionName);
+            regionDirs = regionNames;
+        }
+        else
+        {
+            regionNames = wordList(1, fvMesh::defaultRegion);
+            regionDirs = wordList(1, word::null);
+        }
+    }
+
+    labelList cellToProc;
+
+    forAll(regionNames, regioni)
+    {
+        const word& regionName = regionNames[regioni];
+        const word& regionDir = regionDirs[regioni];
+
+        Info<< "\n\nDecomposing mesh " << regionName << nl << endl;
+        Info<< "Create mesh..." << flush;
+
+        fvMesh mesh
+        (
+            IOobject
+            (
+                regionName,
+                runTime.timeName(),
+                runTime,
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE,
+                false
+            )
+        );
+
+        Info<< " nCells = " << mesh.nCells() << endl;
+
+        // Expected format is a simple ASCII list
+        cellToProc.setSize(mesh.nCells());
+        {
+            IFstream is(decompFile);
+
+            forAll(cellToProc, celli)
+            {
+                cellToProc[celli] = readLabel(is);
+            }
+        }
+
+        const label nDomains = max(cellToProc) + 1;
+
+        CompactListList<label> cellCells;
+        decompositionMethod::calcCellCells
+        (
+            mesh,
+            identity(mesh.nCells()),
+            mesh.nCells(),
+            false,
+            cellCells
+        );
+
+        decompositionInformation info
+        (
+            cellCells,
+            cellToProc,
+            nDomains
+        );
+
+        if (verbose)
+        {
+            info.printDetails(Info);
+            Info<< nl;
+        }
+        info.printSummary(Info);
+    }
+
+    Info<< "\nEnd\n" << endl;
+
+    return 0;
+}
+
+
+// ************************************************************************* //

applications/test/decomposePar/Make/files

@@ -0,0 +1,3 @@
+Test-decomposePar.C
+
+EXE = $(FOAM_APPBIN)/Test-decomposePar

applications/test/decomposePar/Make/options

@@ -0,0 +1,12 @@
+EXE_INC = \
+    -I$(LIB_SRC)/parallel/decompose/decompose/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/regionModels/regionModel/lnInclude
+
+EXE_LIBS = \
+    -ldecompose \
+    -ldecompositionMethods \
+    -L$(FOAM_LIBBIN)/dummy -lmetisDecomp -lscotchDecomp -lkahipDecomp \
+    -lregionModels

applications/test/decomposePar/Test-decomposePar.C

@@ -0,0 +1,278 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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/>.
+
+Application
+    decomposePar
+
+Group
+    grpParallelUtilities
+
+Description
+    Automatically decomposes a mesh and fields of a case for parallel
+    execution of OpenFOAM.
+
+Usage
+    \b decomposePar [OPTION]
+
+    Options:
+      - \par -region \<regionName\>
+        Decompose named region. Does not check for existence of processor*.
+
+      - \par -allRegions
+        Decompose all regions in regionProperties. Does not check for
+        existence of processor*.
+
+      - \par -constant
+
+      - \par -time xxx:yyy
+        Override controlDict settings and decompose selected times. Does not
+        re-decompose the mesh i.e. does not handle moving mesh or changing
+        mesh cases.
+
+\*---------------------------------------------------------------------------*/
+
+#include "OSspecific.H"
+#include "fvCFD.H"
+#include "cpuTime.H"
+#include "IOobjectList.H"
+#include "regionProperties.H"
+#include "decompositionInformation.H"
+#include "decompositionModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+int main(int argc, char *argv[])
+{
+    argList::addNote
+    (
+        "decompose a mesh and fields of a case for parallel execution"
+    );
+
+    argList::noParallel();
+    argList::addOption
+    (
+        "decomposeParDict",
+        "file",
+        "read decomposePar dictionary from specified location"
+    );
+    #include "addRegionOption.H"
+    argList::addBoolOption
+    (
+        "allRegions",
+        "operate on all regions in regionProperties"
+    );
+    argList::addBoolOption
+    (
+        "verbose",
+        "more information about decomposition"
+    );
+    argList::addOption
+    (
+        "domains",
+        "N"
+        "override numberOfSubdomains"
+    );
+
+    argList::addOption
+    (
+        "method",
+        "name"
+        "override method"
+    );
+
+    // Include explicit constant options, have zero from time range
+    timeSelector::addOptions(true, false);
+
+    #include "setRootCase.H"
+
+    const bool region     = args.optionFound("region");
+    const bool allRegions = args.optionFound("allRegions");
+    const bool verbose    = args.optionFound("verbose");
+
+    const label numSubdomains =
+        args.optionLookupOrDefault<label>("domains", 0);
+
+    const word methodName =
+        args.optionLookupOrDefault<word>("method", word::null);
+
+
+    // Set time from database
+    #include "createTime.H"
+    // Allow override of time
+    instantList times = timeSelector::selectIfPresent(runTime, args);
+
+    // Allow override of decomposeParDict location
+    fileName decompDictFile;
+    args.optionReadIfPresent("decomposeParDict", decompDictFile);
+
+    wordList regionNames;
+    wordList regionDirs;
+    if (allRegions)
+    {
+        Info<< "Decomposing all regions in regionProperties" << nl << endl;
+        regionProperties rp(runTime);
+        forAllConstIters(rp, iter)
+        {
+            const wordList& regions = iter();
+            forAll(regions, i)
+            {
+                if (findIndex(regionNames, regions[i]) == -1)
+                {
+                    regionNames.append(regions[i]);
+                }
+            }
+        }
+        regionDirs = regionNames;
+    }
+    else
+    {
+        word regionName;
+        if (args.optionReadIfPresent("region", regionName))
+        {
+            regionNames = wordList(1, regionName);
+            regionDirs = regionNames;
+        }
+        else
+        {
+            regionNames = wordList(1, fvMesh::defaultRegion);
+            regionDirs = wordList(1, word::null);
+        }
+    }
+
+
+    forAll(regionNames, regioni)
+    {
+        const word& regionName = regionNames[regioni];
+        const word& regionDir = regionDirs[regioni];
+
+        Info<< "\n\nDecomposing mesh " << regionName << nl << endl;
+        Info<< "Create mesh..." << flush;
+
+        fvMesh mesh
+        (
+            IOobject
+            (
+                regionName,
+                runTime.timeName(),
+                runTime,
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE,
+                false
+            )
+        );
+
+        Info<< " nCells = " << mesh.nCells() << endl;
+
+        Info<< "\nCalculating distribution of cells" << endl;
+        cpuTime decompositionTime;
+
+        const decompositionModel& model = decompositionModel::New
+        (
+            mesh,
+            decompDictFile
+        );
+
+        // Allow command-line override for quick testing
+
+        dictionary& modelDict = const_cast<decompositionModel&>(model);
+
+        if (numSubdomains)
+        {
+            modelDict.add
+            (
+                word("numberOfSubdomains"),
+                numSubdomains,
+                true
+            );
+        }
+
+        if (!methodName.empty())
+        {
+            modelDict.add
+            (
+                word("method"),
+                methodName,
+                true
+            );
+        }
+
+        scalarField cellWeights;
+        word weightName;
+        if (model.readIfPresent("weightField", weightName))
+        {
+            volScalarField weights
+            (
+                IOobject
+                (
+                    weightName,
+                    mesh.time().timeName(),
+                    mesh,
+                    IOobject::MUST_READ,
+                    IOobject::NO_WRITE
+                ),
+                mesh
+            );
+            cellWeights = weights.primitiveField();
+        }
+
+        decompositionMethod& method = model.decomposer();
+
+        CompactListList<label> cellCells;
+        decompositionMethod::calcCellCells
+        (
+            mesh,
+            identity(mesh.nCells()),
+            mesh.nCells(),
+            false,
+            cellCells
+        );
+
+        labelList cellToProc = method.decompose(mesh, cellWeights);
+
+        Info<< "\nFinished decomposition into "
+            << method.nDomains() << " domains in "
+            << decompositionTime.elapsedCpuTime()
+            << " s" << nl << endl;
+
+        decompositionInformation info
+        (
+            cellCells,
+            cellToProc,
+            method.nDomains()
+        );
+
+        if (verbose)
+        {
+            info.printDetails(Info);
+            Info<< nl;
+        }
+        info.printSummary(Info);
+    }
+
+    Info<< "\nEnd\n" << endl;
+
+    return 0;
+}
+
+
+// ************************************************************************* //

applications/test/foamToMetisGraph/Make/files

@@ -0,0 +1,3 @@
+foamToMetisGraph.C
+
+EXE = $(FOAM_USER_APPBIN)/foamToMetisGraph

applications/test/foamToMetisGraph/Make/options

@@ -0,0 +1,3 @@
+EXE_INC =
+
+EXE_LIBS =

applications/test/foamToMetisGraph/foamToMetisGraph.C

@@ -0,0 +1,86 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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/>.
+
+Description
+    Create a metis graph file representation of an OpenFOAM mesh
+
+\*---------------------------------------------------------------------------*/
+
+#include "argList.H"
+#include "Time.H"
+#include "polyMesh.H"
+#include "OFstream.H"
+
+using namespace Foam;
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// Main program:
+
+int main(int argc, char *argv[])
+{
+    argList::noParallel();
+    argList::noFunctionObjects();
+    argList::addNote
+    (
+        "Create a metis graph file representation for an OpenFOAM mesh"
+    );
+
+    #include "setRootCase.H"
+    #include "createTime.H"
+    #include "createPolyMesh.H"
+
+    const labelListList& cellCells = mesh.cellCells();
+
+    // No. of Nodes = nCells
+    // No. of Edges connecting Nodes = nInternalFaces
+
+    OFstream os(args.caseName() + ".graph", IOstream::ASCII);
+
+    os  << "%% metis graph file, of an OpenFOAM mesh %%" << nl
+        << "%% nCells=" << mesh.nCells()
+        << " nFaces=" << mesh.nFaces()
+        << " nInternalFaces=" << mesh.nInternalFaces() << nl;
+
+    os  << cellCells.size() << " " << mesh.nInternalFaces() << nl;
+
+    for (const auto& edges : cellCells)
+    {
+        forAll(edges, i)
+        {
+            if (i) os << " ";
+            os << edges[i] + 1;  // index starts at 1.
+        }
+        os << nl;
+    }
+
+    Info<<"Wrote graph with "
+        << mesh.nCells() << " nodes and "
+        << mesh.nInternalFaces() << " edges to "
+        << os.name() << nl;
+
+    Info<< nl << "End\n" << endl;
+
+    return 0;
+}
+
+// ************************************************************************* //

applications/utilities/mesh/manipulation/renumberMesh/Make/options

@@ -18,4 +18,4 @@ EXE_LIBS = \
     -lreconstruct \
     $(LINK_FLAGS) \
     -ldecompositionMethods \
-    -L$(FOAM_LIBBIN)/dummy -lmetisDecomp -lscotchDecomp
+    -L$(FOAM_LIBBIN)/dummy -lmetisDecomp -lscotchDecomp -lkahipDecomp

applications/utilities/parallelProcessing/decomposePar/Make/options

@@ -11,7 +11,8 @@ EXE_LIBS = \
     -ldynamicMesh \
     -ldecompose \
     -lgenericPatchFields \
-    -ldecompositionMethods -L$(FOAM_LIBBIN)/dummy -lmetisDecomp -lscotchDecomp \
+    -ldecompositionMethods \
+    -L$(FOAM_LIBBIN)/dummy -lmetisDecomp -lscotchDecomp -lkahipDecomp \
     -llagrangian \
     -ldynamicMesh \
     -lregionModels

etc/config.sh/kahip

@@ -0,0 +1,59 @@
+#----------------------------------*-sh-*--------------------------------------
+# =========                 |
+# \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+#  \\    /   O peration     |
+#   \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+#    \\/     M anipulation  |
+#------------------------------------------------------------------------------
+# This file is part of OpenFOAM, licensed under the GNU General Public License
+# <http://www.gnu.org/licenses/>.
+#
+# File
+#     etc/config.sh/kahip
+#
+# Description
+#     Setup for kahip include/libraries.
+#     Sourced during wmake process only.
+#
+#     Normally used to specify the kahip version and location for a
+#     ThirdParty installation.
+#
+#     If using system-wide installations, use the following setting:
+#
+#         KAHIP_VERSION=kahip-system
+#
+#     If the system kahip is unusable (eg, too old) and you don't
+#     have or want a ThirdParty installation:
+#
+#         KAHIP_VERSION=kahip-none
+#
+#     If using a central installation, but not located under ThirdParty:
+#         - specify kahip-system
+#         - provide full path for KAHIP_ARCH_PATH
+#
+# Note
+#     A csh version is not needed, since the values here are only sourced
+#     during the wmake process.
+#
+#     KaHIP can also be entirely disabled, by either renaming this file or
+#     by creating an empty one with the same name at a user or site location.
+#
+#     KaHIP is 32-bit precision only.
+#     An Int64 OpenFOAM version can use it, but the mesh size is limited
+#     accordingly.
+#
+#     If KaHIP was compiled with openmp, you may need to add in additional
+#     compile or link flags in KAHIP_COMP_FLAGS KAHIP_LINK_FLAGS
+#
+#------------------------------------------------------------------------------
+# USER EDITABLE PART: Changes made here may be lost with the next upgrade
+
+KAHIP_VERSION=kahip-2.00
+export KAHIP_ARCH_PATH=$WM_THIRD_PARTY_DIR/platforms/$WM_ARCH$WM_COMPILER/$KAHIP_VERSION
+
+# Adjust as required
+# export KAHIP_COMP_FLAGS="-fopenmp"
+# export KAHIP_LINK_FLAGS="-lgomp"
+
+# END OF (NORMAL) USER EDITABLE PART
+#------------------------------------------------------------------------------

etc/config.sh/unset

@@ -134,6 +134,7 @@ unset GPERFTOOLS_ARCH_PATH
 unset GMP_ARCH_PATH
 unset MPFR_ARCH_PATH
 unset MESA_ARCH_PATH
+unset KAHIP_ARCH_PATH
 unset METIS_ARCH_PATH
 unset SCOTCH_ARCH_PATH
 

src/dummyThirdParty/Allwmake

@@ -7,6 +7,7 @@ cd ${0%/*} || exit 1    # Run from this directory
 wmake $targetType scotchDecomp
 wmake $targetType ptscotchDecomp
 wmake $targetType metisDecomp
+wmake $targetType kahipDecomp
 wmake $targetType MGridGen
 
 #------------------------------------------------------------------------------

src/dummyThirdParty/kahipDecomp/Make/files

@@ -0,0 +1,3 @@
+dummyKahipDecomp.C
+
+LIB = $(FOAM_LIBBIN)/dummy/libkahipDecomp

src/dummyThirdParty/kahipDecomp/Make/options

@@ -0,0 +1,5 @@
+EXE_INC = \
+    -I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
+    -I$(LIB_SRC)/parallel/decompose/kahipDecomp/lnInclude
+
+LIB_LIBS =

src/dummyThirdParty/kahipDecomp/dummyKahipDecomp.C

@@ -0,0 +1,83 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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 "kahipDecomp.H"
+#include "addToRunTimeSelectionTable.H"
+#include "Time.H"
+
+static const char* notImplementedMessage =
+"You are trying to use kahip but do not have the kahipDecomp library loaded."
+"\nThis message is from the dummy kahipDecomp stub library instead.\n"
+"\n"
+"Please install kahip and make sure that libkahip.so is in your "
+"LD_LIBRARY_PATH.\n"
+"The kahipDecomp library can then be built from "
+"src/parallel/decompose/kahipDecomp and dynamically loading or linking"
+" this library will add kahip as a decomposition method.\n";
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(kahipDecomp, 0);
+
+    addToRunTimeSelectionTable
+    (
+        decompositionMethod,
+        kahipDecomp,
+        dictionary
+    );
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+Foam::label Foam::kahipDecomp::decomposeSerial
+(
+    const UList<label>& adjncy,
+    const UList<label>& xadj,
+    const UList<scalar>& cellWeights,
+    List<label>& decomp
+)
+{
+    FatalErrorInFunction
+        << notImplementedMessage << exit(FatalError);
+
+    return -1;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::kahipDecomp::kahipDecomp
+(
+    const dictionary& decompositionDict
+)
+:
+    metisLikeDecomp(decompositionDict)
+{}
+
+
+// ************************************************************************* //

src/dummyThirdParty/kahipDecomp/kaHIP_interface.h

@@ -0,0 +1,37 @@
+#ifndef KAHIP_H
+#define KAHIP_H
+
+/* *** DUMMY VERSION of kaHIP_interface.h - this file should not be included if you have KaHIP
+ *     installed in the correct position in $WM_THIRD_PARTY_DIR - see
+ *     decompositionMethods/kahipDecomp/Make/options
+ */
+
+#warning "Dummy kahip.h - included since it cannot find KaHIP installation."
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+void kaffpa
+(
+    // [inputs]
+    int* n,
+    int* vwgt,
+    int* xadj,
+    int* adjcwgt,
+    int* adjncy,
+    int* nparts,
+    double* imbalance,
+    bool suppress_output,
+    int seed,
+    int mode,
+    // [outputs]
+    int* edgecut,
+    int* part
+);
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif

src/dummyThirdParty/metisDecomp/dummyMetisDecomp.C

@@ -35,8 +35,7 @@ static const char* notImplementedMessage =
 "LD_LIBRARY_PATH.\n"
 "The metisDecomp library can then be built from "
 "src/parallel/decompose/metisDecomp and dynamically loading or linking"
-" this library will add metis as a decomposition method.\n"
+" this library will add metis as a decomposition method.\n";
-"Please be aware that there are license restrictions on using Metis.";
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -55,12 +54,12 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-Foam::label Foam::metisDecomp::decompose
+Foam::label Foam::metisDecomp::decomposeSerial
 (
-    const List<label>& adjncy,
+    const UList<label>& adjncy,
-    const List<label>& xadj,
+    const UList<label>& xadj,
-    const scalarField& cellWeights,
+    const UList<scalar>& cellWeights,
-    List<label>& finalDecomp
+    List<label>& decomp
 )
 {
     FatalErrorInFunction
@@ -77,53 +76,8 @@ Foam::metisDecomp::metisDecomp
     const dictionary& decompositionDict
 )
 :
-    decompositionMethod(decompositionDict)
+    metisLikeDecomp(decompositionDict)
 {}
 
 
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::labelList Foam::metisDecomp::decompose
-(
-    const polyMesh& mesh,
-    const pointField& points,
-    const scalarField& pointWeights
-)
-{
-    FatalErrorInFunction
-        << notImplementedMessage << exit(FatalError);
-
-    return labelList();
-}
-
-
-Foam::labelList Foam::metisDecomp::decompose
-(
-    const polyMesh& mesh,
-    const labelList& agglom,
-    const pointField& agglomPoints,
-    const scalarField& agglomWeights
-)
-{
-    FatalErrorInFunction
-        << notImplementedMessage << exit(FatalError);
-
-    return labelList();
-}
-
-
-Foam::labelList Foam::metisDecomp::decompose
-(
-    const labelListList& globalCellCells,
-    const pointField& cellCentres,
-    const scalarField& cellWeights
-)
-{
-    FatalErrorInFunction
-        << notImplementedMessage << exit(FatalError);
-
-    return labelList();
-}
-
-
 // ************************************************************************* //

src/dummyThirdParty/ptscotchDecomp/dummyPtscotchDecomp.C

@@ -54,17 +54,19 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+void Foam::ptscotchDecomp::graphPath(const polyMesh& unused)
+{}
+
+
 void Foam::ptscotchDecomp::check(const int retVal, const char* str)
 {}
 
 
 Foam::label Foam::ptscotchDecomp::decompose
 (
-    const fileName& meshPath,
+    const UList<label>& initxadj,
-    const List<label>& initxadj,
+    const UList<label>& initadjncy,
-    const List<label>& initadjncy,
+    const UList<scalar>& initcWeights,
-    const scalarField& initcWeights,
-
     List<label>& finalDecomp
 ) const
 {
@@ -77,12 +79,11 @@ Foam::label Foam::ptscotchDecomp::decompose
 
 Foam::label Foam::ptscotchDecomp::decompose
 (
-    const fileName& meshPath,
     const label adjncySize,
     const label adjncy[],
     const label xadjSize,
     const label xadj[],
-    const scalarField& cWeights,
+    const UList<scalar>& cWeights,
     List<label>& finalDecomp
 ) const
 {

src/dummyThirdParty/scotchDecomp/dummyScotchDecomp.C

@@ -53,18 +53,20 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+void Foam::scotchDecomp::graphPath(const polyMesh& unused)
+{}
+
+
 void Foam::scotchDecomp::check(const int retVal, const char* str)
 {}
 
 
-Foam::label Foam::scotchDecomp::decompose
+Foam::label Foam::scotchDecomp::decomposeSerial
 (
-    const fileName& meshPath,
+    const UList<label>& adjncy,
-    const List<label>& adjncy,
+    const UList<label>& xadj,
-    const List<label>& xadj,
+    const UList<scalar>& cWeights,
-    const scalarField& cWeights,
+    List<label>& decomp
-
-    List<label>& finalDecomp
 )
 {
     FatalErrorInFunction
@@ -81,7 +83,7 @@ Foam::scotchDecomp::scotchDecomp
     const dictionary& decompositionDict
 )
 :
-    decompositionMethod(decompositionDict)
+    metisLikeDecomp(decompositionDict)
 {}
 
 

src/parallel/decompose/Allwclean

@@ -17,8 +17,9 @@ wcleanMpiLib()
 }
 
 
-wclean scotchDecomp
 wclean metisDecomp
+wclean kahipDecomp
+wclean scotchDecomp
 wclean decompositionMethods
 wclean decompose
 wcleanMpiLib ptscotchDecomp

src/parallel/decompose/Allwmake

@@ -25,6 +25,54 @@ findFirstFile()
 }
 
 
+# Test for kahip.
+# - return 0 and export KAHIP_ARCH_PATH on success
+hasKahip()
+{
+    local warning="==> skip kahip"
+    local header label settings
+
+    unset KAHIP_ARCH_PATH KAHIP_VERSION
+    unset KAHIP_COMP_FLAGS KAHIP_LINK_FLAGS
+    settings=$($WM_PROJECT_DIR/bin/foamEtcFile config.sh/kahip) || {
+        echo "$warning (no config.sh/kahip settings)"
+        return 1
+    }
+
+    . $settings
+    if [ -z "$KAHIP_ARCH_PATH" -o "${KAHIP_ARCH_PATH##*-}" = none ]
+    then
+        echo "$warning (not available)"
+        return 1
+    fi
+
+    # Header
+    header=$(findFirstFile \
+        $KAHIP_ARCH_PATH/include/kaHIP_interface.h \
+        /usr/include/kaHIP_interface.h \
+    ) || {
+        echo "$warning (no header)"
+        return 2  # file not found
+    }
+
+    # Library
+    [ "${KAHIP_ARCH_PATH##*-}" = system ] || \
+    findFirstFile \
+        $KAHIP_ARCH_PATH/lib/libkahip.a  \
+        $KAHIP_ARCH_PATH/lib$WM_COMPILER_LIB_ARCH/libkahip.a \
+    > /dev/null || {
+        echo "$warning (missing library)"
+        return 2
+    }
+
+    # kahip itself is 32-bit int, but our interface itself handles some
+    # 64-bit conversion (mesh size).
+
+    export KAHIP_ARCH_PATH
+    echo "kahip (label=32) - $KAHIP_ARCH_PATH"
+}
+
+
 # Test for metis.
 # - return 0 and export METIS_ARCH_PATH on success
 hasMetis()
@@ -208,6 +256,11 @@ then
     wmake $targetType metisDecomp
 fi
 
+if hasKahip
+then
+    wmake $targetType kahipDecomp
+fi
+
 wmake $targetType decompositionMethods
 wmake $targetType decompose
 

src/parallel/decompose/decompose/Make/files

@@ -1,3 +1,4 @@
+decompositionInformation.C
 decompositionModel.C
 fvFieldDecomposer.C
 

src/parallel/decompose/decompose/decompositionInformation.C

@@ -0,0 +1,236 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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 "decompositionInformation.H"
+#include "ListOps.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+void Foam::decompositionInformation::populate
+(
+    const labelUList& adjncy,
+    const labelUList& xadj,
+    const labelUList& decomp,
+    const label nDomain
+)
+{
+    nDomains_ = nDomain;
+
+    distrib_.clear();
+    distrib_.setSize(nDomain);
+
+    for (labelList& subdist : distrib_)
+    {
+        subdist.clear();
+        subdist.setSize(nDomain, Zero);
+    }
+
+    const label nCells = xadj.size()-1;
+    for (label celli = 0; celli < nCells; ++celli)
+    {
+        const label ownProc = decomp[celli];
+
+        labelList& subdist = distrib_[ownProc];
+
+        // Number of cells
+        ++subdist[ownProc];
+
+        for (label i = xadj[celli]; i < xadj[celli+1]; ++i)
+        {
+            const label neiProc = decomp[adjncy[i]];
+
+            if (neiProc != ownProc)
+            {
+                // Number of processor faces
+                ++subdist[neiProc];
+            }
+        }
+    }
+
+    // Build summary
+
+    labelList cellsCount(nDomains_, Zero);
+    labelList neighCount(nDomains_, Zero);
+    labelList facesCount(nDomains_, Zero);
+
+    forAll(distrib_, ownProc)
+    {
+        const labelList& subdist = distrib_[ownProc];
+
+        cellsCount[ownProc] = subdist[ownProc];
+
+        forAll(subdist, neiProc)
+        {
+            const label n = subdist[neiProc];
+
+            if (n && ownProc != neiProc)
+            {
+                ++neighCount[ownProc];
+                facesCount[ownProc] += n;
+            }
+        }
+    }
+
+    const label n2 = (nDomains_ / 2);
+
+    sort(cellsCount);
+    cellsInfo_.min    = cellsCount.first();
+    cellsInfo_.max    = cellsCount.last();
+    cellsInfo_.median = cellsCount[n2];
+
+    sort(neighCount);
+    neighInfo_.min    = neighCount.first();
+    neighInfo_.max    = neighCount.last();
+    neighInfo_.median = neighCount[n2];
+
+    sort(facesCount);
+    facesInfo_.min    = facesCount.first();
+    facesInfo_.max    = facesCount.last();
+    facesInfo_.median = facesCount[n2];
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::decompositionInformation::decompositionInformation
+(
+    const labelUList& adjncy,
+    const labelUList& xadj,
+    const labelUList& decomp,
+    const label nDomains
+)
+:
+    distrib_(),
+    nDomains_(0)
+{
+    populate(adjncy, xadj, decomp, nDomains);
+}
+
+
+Foam::decompositionInformation::decompositionInformation
+(
+    const CompactListList<label>& cellCells,
+    const labelUList& decomp,
+    const label nDomains
+)
+:
+    distrib_(),
+    nDomains_(0)
+{
+    populate(cellCells.m(), cellCells.offsets(), decomp, nDomains);
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+void Foam::decompositionInformation::clear()
+{
+    distrib_.clear();
+    cellsInfo_.clear();
+    neighInfo_.clear();
+    facesInfo_.clear();
+}
+
+
+void Foam::decompositionInformation::printSummary(Ostream& os) const
+{
+    os  << "Cells "; cellsInfo_.print(os) << nl;
+    os  << "Neigh "; neighInfo_.print(os)<< nl;
+    os  << "Faces "; facesInfo_.print(os)<< nl;
+}
+
+
+void Foam::decompositionInformation::printDetails(Ostream& os) const
+{
+    forAll(distrib_, ownProc)
+    {
+        const labelList& subdist = distrib_[ownProc];
+
+        // First pass:
+        label neighCount = 0;
+        label facesCount = 0;
+
+        forAll(subdist, neiProc)
+        {
+            const label n = subdist[neiProc];
+
+            if (n && ownProc != neiProc)
+            {
+                ++neighCount;
+                facesCount += n;
+            }
+        }
+
+        os  << "Part[" << ownProc << "] cells:" << subdist[ownProc]
+            << " neigh:" << neighCount
+            << " faces:" << facesCount;
+
+        // Second pass with details:
+        if (facesCount)
+        {
+            Info<< " ";
+
+            forAll(subdist, neiProc)
+            {
+                const label n = subdist[neiProc];
+
+                if (n && ownProc != neiProc)
+                {
+                    os  << " (" << neiProc << " " << n << ")";
+                }
+            }
+        }
+
+        os  << nl;
+    }
+}
+
+
+void Foam::decompositionInformation::printAll(Ostream& os) const
+{
+    printDetails(os);
+    printSummary(os);
+}
+
+
+// * * * * * * * * * * * * * * * Ostream Operator  * * * * * * * * * * * * * //
+
+
+Foam::Ostream& Foam::decompositionInformation::stats::print(Ostream& os) const
+{
+    os  << "max/median/min: "
+        << this->max << " / " << this->median << " / " << this->min;
+
+    if (this->median)
+    {
+        const scalar ratio = scalar(100*this->max)/this->median;
+
+        os  << " (" << ratio << "%)";
+    }
+
+    return os;
+}
+
+// ************************************************************************* //

src/parallel/decompose/decompose/decompositionInformation.H

@@ -0,0 +1,149 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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::decompositionInformation
+
+Description
+    Abstract base class for decomposition
+
+SourceFiles
+    decompositionInformation.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef decompositionInformation_H
+#define decompositionInformation_H
+
+#include "polyMesh.H"
+#include "CompactListList.H"
+#include "IOstreams.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                  Class decompositionInformation Declaration
+\*---------------------------------------------------------------------------*/
+
+class decompositionInformation
+{
+    //- Simple storage for organizing min/max/median quantities
+    struct stats
+    {
+        label min;
+        label max;
+        label median;
+
+        stats() : min(0), max(0), median(0) {}
+
+        void clear()
+        {
+            min = 0;
+            max = 0;
+            median = 0;
+        }
+
+        Ostream& print(Ostream& os) const;
+    };
+
+    // Private data
+
+        labelListList distrib_;
+        label nDomains_;
+
+        stats cellsInfo_;
+        stats neighInfo_;
+        stats facesInfo_;
+
+
+    // Private Member Functions
+
+        //- Populate from cell decomposition list
+        void populate
+        (
+            const labelUList& adjncy,
+            const labelUList& xadj,
+            const labelUList& decomp,
+            const label nDomains
+        );
+
+        //- Disallow default bitwise copy construct and assignment
+        decompositionInformation(const decompositionInformation&) = delete;
+        void operator=(const decompositionInformation&) = delete;
+
+
+public:
+
+    // Constructors
+
+        //- Construct
+        decompositionInformation
+        (
+            const labelUList& adjncy,
+            const labelUList& xadj,
+            const labelUList& decomp,
+            const label nDomains
+        );
+
+        //- Construct
+        decompositionInformation
+        (
+            const CompactListList<label>& cellCells,
+            const labelUList& decomp,
+            const label nDomains
+        );
+
+
+    //- Destructor
+    ~decompositionInformation()
+    {}
+
+
+    // Member Functions
+
+        void clear();
+
+        label nDomains() const
+        {
+            return nDomains_;
+        }
+
+
+        void printSummary(Ostream& os) const;
+
+        void printDetails(Ostream& os) const;
+
+        void printAll(Ostream& os) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/parallel/decompose/decompositionMethods/Make/files

@@ -4,6 +4,7 @@ simpleGeomDecomp/simpleGeomDecomp.C
 hierarchGeomDecomp/hierarchGeomDecomp.C
 manualDecomp/manualDecomp.C
 multiLevelDecomp/multiLevelDecomp.C
+metisLikeDecomp/metisLikeDecomp.C
 structuredDecomp/structuredDecomp.C
 noDecomp/noDecomp.C
 

src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C

@@ -21,9 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-InClass
-    decompositionMethod
-
 \*---------------------------------------------------------------------------*/
 
 #include "decompositionMethod.H"
@@ -49,6 +46,7 @@ namespace Foam
     defineRunTimeSelectionTable(decompositionMethod, dictionary);
 }
 
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::decompositionMethod::decompositionMethod
@@ -349,8 +347,8 @@ void Foam::decompositionMethod::calcCellCells
 
     for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
     {
-        label own = agglom[faceOwner[facei]];
+        const label own = agglom[faceOwner[facei]];
-        label nei = agglom[faceNeighbour[facei]];
+        const label nei = agglom[faceNeighbour[facei]];
 
         nFacesPerCell[own]++;
         nFacesPerCell[nei]++;
@@ -367,7 +365,7 @@ void Foam::decompositionMethod::calcCellCells
 
             forAll(pp, i)
             {
-                label own = agglom[faceOwner[facei]];
+                const label own = agglom[faceOwner[facei]];
 
                 label globalNei = globalNeighbour[bFacei];
                 if
@@ -399,8 +397,8 @@ void Foam::decompositionMethod::calcCellCells
     // For internal faces is just offsetted owner and neighbour
     for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
     {
-        label own = agglom[faceOwner[facei]];
+        const label own = agglom[faceOwner[facei]];
-        label nei = agglom[faceNeighbour[facei]];
+        const label nei = agglom[faceNeighbour[facei]];
 
         m[offsets[own] + nFacesPerCell[own]++] = globalAgglom.toGlobal(nei);
         m[offsets[nei] + nFacesPerCell[nei]++] = globalAgglom.toGlobal(own);
@@ -418,9 +416,9 @@ void Foam::decompositionMethod::calcCellCells
 
             forAll(pp, i)
             {
-                label own = agglom[faceOwner[facei]];
+                const label own = agglom[faceOwner[facei]];
 
-                label globalNei = globalNeighbour[bFacei];
+                const label globalNei = globalNeighbour[bFacei];
 
                 if
                 (
@@ -457,7 +455,7 @@ void Foam::decompositionMethod::calcCellCells
         nbrCells.clear();
         nbrCells.insert(globalAgglom.toGlobal(celli));
 
-        label endIndex = cellCells.offsets()[celli+1];
+        const label endIndex = cellCells.offsets()[celli+1];
 
         for (label i = startIndex; i < endIndex; i++)
         {
@@ -555,8 +553,8 @@ void Foam::decompositionMethod::calcCellCells
 
     for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
     {
-        label own = agglom[faceOwner[faceI]];
+        const label own = agglom[faceOwner[faceI]];
-        label nei = agglom[faceNeighbour[faceI]];
+        const label nei = agglom[faceNeighbour[faceI]];
 
         nFacesPerCell[own]++;
         nFacesPerCell[nei]++;
@@ -573,9 +571,9 @@ void Foam::decompositionMethod::calcCellCells
 
             forAll(pp, i)
             {
-                label own = agglom[faceOwner[faceI]];
+                const label own = agglom[faceOwner[faceI]];
 
-                label globalNei = globalNeighbour[bFaceI];
+                const label globalNei = globalNeighbour[bFaceI];
                 if
                 (
                    !globalAgglom.isLocal(globalNei)
@@ -607,11 +605,11 @@ void Foam::decompositionMethod::calcCellCells
     // For internal faces is just offsetted owner and neighbour
     for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
     {
-        label own = agglom[faceOwner[faceI]];
+        const label own = agglom[faceOwner[faceI]];
-        label nei = agglom[faceNeighbour[faceI]];
+        const label nei = agglom[faceNeighbour[faceI]];
 
-        label ownIndex = offsets[own] + nFacesPerCell[own]++;
+        const label ownIndex = offsets[own] + nFacesPerCell[own]++;
-        label neiIndex = offsets[nei] + nFacesPerCell[nei]++;
+        const label neiIndex = offsets[nei] + nFacesPerCell[nei]++;
 
         m[ownIndex] = globalAgglom.toGlobal(nei);
         w[ownIndex] = mag(mesh.faceAreas()[faceI]);
@@ -631,9 +629,9 @@ void Foam::decompositionMethod::calcCellCells
 
             forAll(pp, i)
             {
-                label own = agglom[faceOwner[faceI]];
+                const label own = agglom[faceOwner[faceI]];
 
-                label globalNei = globalNeighbour[bFaceI];
+                const label globalNei = globalNeighbour[bFaceI];
 
                 if
                 (
@@ -672,7 +670,7 @@ void Foam::decompositionMethod::calcCellCells
         nbrCells.clear();
         nbrCells.insert(globalAgglom.toGlobal(cellI));
 
-        label endIndex = cellCells.offsets()[cellI+1];
+        const label endIndex = cellCells.offsets()[cellI+1];
 
         for (label i = startIndex; i < endIndex; i++)
         {

src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.H

@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2015-2017 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.

src/parallel/decompose/decompositionMethods/metisLikeDecomp/metisLikeDecomp.C

@@ -0,0 +1,264 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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 "metisLikeDecomp.H"
+#include "Time.H"
+#include "globalIndex.H"
+
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
+
+Foam::label Foam::metisLikeDecomp::decomposeGeneral
+(
+    const UList<label>& adjncy,
+    const UList<label>& xadj,
+    const UList<scalar>& cWeights,
+    List<label>& decomp
+)
+{
+    if (!Pstream::parRun())
+    {
+        return decomposeSerial
+        (
+            adjncy,
+            xadj,
+            cWeights,
+            decomp
+        );
+    }
+
+    if (debug)
+    {
+        Info<< type() << "Decomp : running in parallel."
+            << " Decomposing all of graph on master processor." << endl;
+    }
+    globalIndex globalCells(xadj.size()-1);
+    label nTotalConnections = returnReduce(adjncy.size(), sumOp<label>());
+
+    // Send all to master. Use scheduled to save some storage.
+    if (Pstream::master())
+    {
+        List<label> allAdjncy(nTotalConnections);
+        List<label> allXadj(globalCells.size()+1);
+        List<scalar> allWeights(globalCells.size());
+
+        // Insert my own
+        label nTotalCells = 0;
+        forAll(cWeights, celli)
+        {
+            allXadj[nTotalCells] = xadj[celli];
+            allWeights[nTotalCells++] = cWeights[celli];
+        }
+        nTotalConnections = 0;
+        forAll(adjncy, i)
+        {
+            allAdjncy[nTotalConnections++] = adjncy[i];
+        }
+
+        for (int slave=1; slave<Pstream::nProcs(); ++slave)
+        {
+            IPstream fromSlave(Pstream::commsTypes::scheduled, slave);
+            List<label> nbrAdjncy(fromSlave);
+            List<label> nbrXadj(fromSlave);
+            List<scalar> nbrWeights(fromSlave);
+
+            // Append.
+            forAll(nbrXadj, celli)
+            {
+                allXadj[nTotalCells] = nTotalConnections+nbrXadj[celli];
+                allWeights[nTotalCells++] = nbrWeights[celli];
+            }
+            // No need to renumber xadj since already global.
+            forAll(nbrAdjncy, i)
+            {
+                allAdjncy[nTotalConnections++] = nbrAdjncy[i];
+            }
+        }
+        allXadj[nTotalCells] = nTotalConnections;
+
+        labelList allDecomp;
+        decomposeSerial
+        (
+            allAdjncy,
+            allXadj,
+            allWeights,
+            allDecomp
+        );
+
+
+        // Send allFinalDecomp back
+        for (int slave=1; slave<Pstream::nProcs(); ++slave)
+        {
+            OPstream toSlave(Pstream::commsTypes::scheduled, slave);
+            toSlave << SubList<label>
+            (
+                allDecomp,
+                globalCells.localSize(slave),
+                globalCells.offset(slave)
+            );
+        }
+
+        // Get my own part (always first)
+        decomp = SubList<label>(allDecomp, globalCells.localSize());
+    }
+    else
+    {
+        // Send my part of the graph (already in global numbering)
+        {
+            OPstream toMaster
+            (
+                Pstream::commsTypes::scheduled,
+                Pstream::masterNo()
+            );
+            toMaster<< adjncy << SubList<label>(xadj, xadj.size()-1)
+                << cWeights;
+        }
+
+        // Receive back decomposition
+        IPstream fromMaster
+        (
+            Pstream::commsTypes::scheduled,
+            Pstream::masterNo()
+        );
+        fromMaster >> decomp;
+    }
+
+    return 0;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::metisLikeDecomp::metisLikeDecomp(const dictionary& decompositionDict)
+:
+    decompositionMethod(decompositionDict)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::labelList Foam::metisLikeDecomp::decompose
+(
+    const polyMesh& mesh,
+    const pointField& points,
+    const scalarField& pointWeights
+)
+{
+    if (points.size() != mesh.nCells())
+    {
+        FatalErrorInFunction
+            << "Can use this decomposition method only for entire mesh" << nl
+            << "and supply one coordinate (cellCentre) for every cell." << nl
+            << "The number of coordinates " << points.size() << nl
+            << "The number of cells in the mesh " << mesh.nCells()
+            << exit(FatalError);
+    }
+
+    CompactListList<label> cellCells;
+    calcCellCells
+    (
+        mesh,
+        identity(mesh.nCells()),
+        mesh.nCells(),
+        false,
+        cellCells
+    );
+
+    // Decompose using default weights
+    labelList decomp;
+    decomposeGeneral(cellCells.m(), cellCells.offsets(), pointWeights, decomp);
+
+    return decomp;
+}
+
+
+Foam::labelList Foam::metisLikeDecomp::decompose
+(
+    const polyMesh& mesh,
+    const labelList& agglom,
+    const pointField& agglomPoints,
+    const scalarField& agglomWeights
+)
+{
+    if (agglom.size() != mesh.nCells())
+    {
+        FatalErrorInFunction
+            << "Size of cell-to-coarse map " << agglom.size()
+            << " differs from number of cells in mesh " << mesh.nCells()
+            << exit(FatalError);
+    }
+
+    // Make Metis CSR (Compressed Storage Format) storage
+    //   adjncy      : contains neighbours (= edges in graph)
+    //   xadj(celli) : start of information in adjncy for celli
+
+    CompactListList<label> cellCells;
+    calcCellCells(mesh, agglom, agglomPoints.size(), false, cellCells);
+
+    // Decompose using default weights
+    labelList decomp;
+    decomposeGeneral(cellCells.m(), cellCells.offsets(), agglomWeights, decomp);
+
+
+    // Rework back into decomposition for original mesh
+    labelList fineDistribution(agglom.size());
+
+    forAll(fineDistribution, i)
+    {
+        fineDistribution[i] = decomp[agglom[i]];
+    }
+
+    return fineDistribution;
+}
+
+
+Foam::labelList Foam::metisLikeDecomp::decompose
+(
+    const labelListList& globalCellCells,
+    const pointField& cellCentres,
+    const scalarField& cellWeights
+)
+{
+    if (cellCentres.size() != globalCellCells.size())
+    {
+        FatalErrorInFunction
+            << "Inconsistent number of cells (" << globalCellCells.size()
+            << ") and number of cell centres (" << cellCentres.size()
+            << ")." << exit(FatalError);
+    }
+
+    // Make Metis CSR (Compressed Storage Format) storage
+    //   adjncy      : contains neighbours (= edges in graph)
+    //   xadj(celli) : start of information in adjncy for celli
+
+    CompactListList<label> cellCells(globalCellCells);
+
+    // Decompose using default weights
+    labelList decomp;
+    decomposeGeneral(cellCells.m(), cellCells.offsets(), cellWeights, decomp);
+
+    return decomp;
+}
+
+// ************************************************************************* //

src/parallel/decompose/decompositionMethods/metisLikeDecomp/metisLikeDecomp.H

@@ -0,0 +1,151 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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::metisLikeDecomp
+
+Description
+    Domain decomposition using METIS-like data structures.
+
+    When run in parallel will collect the entire graph on to the master,
+    decompose and send back.
+
+SourceFiles
+    metisLikeDecomp.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef metisLikeDecomp_H
+#define metisLikeDecomp_H
+
+#include "decompositionMethod.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                       Class metisLikeDecomp Declaration
+\*---------------------------------------------------------------------------*/
+
+class metisLikeDecomp
+:
+    public decompositionMethod
+{
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct and assignment
+        void operator=(const metisLikeDecomp&) = delete;
+        metisLikeDecomp(const metisLikeDecomp&) = delete;
+
+
+protected:
+
+    // Protected Member Functions
+
+        //- Serial and/or collect/distribute for parallel operation
+        virtual label decomposeGeneral
+        (
+            const UList<label>& adjncy,
+            const UList<label>& xadj,
+            const UList<scalar>& cellWeights,
+            List<label>& decomp
+        );
+
+        //- Decomposition with metis-like parameters
+        virtual label decomposeSerial
+        (
+            const UList<label>& adjncy,
+            const UList<label>& xadj,
+            const UList<scalar>& cellWeights,
+            List<label>& decomp
+        ) = 0;
+
+public:
+
+    // Constructors
+
+        //- Construct given the decomposition dictionary
+        metisLikeDecomp(const dictionary& decompositionDict);
+
+
+    //- Destructor
+    virtual ~metisLikeDecomp()
+    {}
+
+
+    // Member Functions
+
+        //- Inherit decompose from decompositionMethod
+        using decompositionMethod::decompose;
+
+        //- Return for every coordinate the wanted processor number.
+        //  Uses the mesh connectivity (if needed).
+        //  Weights get normalised so the minimum value is 1 before truncation
+        //  to an integer so the weights should be multiples of the minimum
+        //  value. The overall sum of weights might otherwise overflow.
+        virtual labelList decompose
+        (
+            const polyMesh& mesh,
+            const pointField& points,
+            const scalarField& pointWeights
+        );
+
+        //- Return for every coordinate the wanted processor number.
+        //  Gets passed agglomeration map (from fine to coarse cells) and coarse
+        //  cell location. Can be overridden by decomposers that provide this
+        //  functionality natively.
+        //  See note on weights above.
+        virtual labelList decompose
+        (
+            const polyMesh& mesh,
+            const labelList& agglom,
+            const pointField& regionPoints,
+            const scalarField& regionWeights
+        );
+
+        //- Return for every coordinate the wanted processor number.
+        //  Explicitly provided mesh connectivity.
+        //  The connectivity is equal to mesh.cellCells() except for
+        //  - in parallel the cell numbers are global cell numbers (starting
+        //    from 0 at processor0 and then incrementing all through the
+        //    processors)
+        //  - the connections are across coupled patches
+        //  See note on weights above.
+        virtual labelList decompose
+        (
+            const labelListList& globalCellCells,
+            const pointField& cellCentres,
+            const scalarField& cellWeights
+        );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/parallel/decompose/kahipDecomp/Make/files

@@ -0,0 +1,3 @@
+kahipDecomp.C
+
+LIB = $(FOAM_LIBBIN)/libkahipDecomp

src/parallel/decompose/kahipDecomp/Make/options

@@ -0,0 +1,13 @@
+EXE_INC = \
+    -I$(KAHIP_ARCH_PATH)/include \
+    -I../decompositionMethods/lnInclude \
+    $(KAHIP_COMP_FLAGS)
+
+/*
+ * The $(KAHIP_ARCH_PATH)/lib$WM_COMPILER_LIB_ARCH path is provided
+ * to support central, non-thirdparty installations
+ */
+LIB_LIBS = \
+    -L$(KAHIP_ARCH_PATH)/lib \
+    -L$(KAHIP_ARCH_PATH)/lib$(WM_COMPILER_LIB_ARCH) \
+    -L$(FOAM_EXT_LIBBIN) $(KAHIP_LINK_FLAGS) -lkahip

src/parallel/decompose/kahipDecomp/kahipDecomp.C

@@ -0,0 +1,329 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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 "kahipDecomp.H"
+#include "addToRunTimeSelectionTable.H"
+#include "Time.H"
+
+#include "kaHIP_interface.h"
+
+#include <string>
+#include <map>
+#include <vector>
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(kahipDecomp, 0);
+    addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary);
+}
+
+
+const Foam::Enum
+<
+    Foam::kahipDecomp::configs
+>
+Foam::kahipDecomp::configNames
+{
+    { kahipDecomp::configs::FAST, "fast" },
+    { kahipDecomp::configs::ECO, "eco" },
+    { kahipDecomp::configs::STRONG, "strong" },
+    { kahipDecomp::configs::FASTSOCIAL, "fast-social" },
+    { kahipDecomp::configs::ECOSOCIAL, "eco-social" },
+    { kahipDecomp::configs::STRONGSOCIAL, "strong-social" },
+};
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+Foam::label Foam::kahipDecomp::decomposeSerial
+(
+    const UList<label>& adjncy,
+    const UList<label>& xadj,
+    const UList<scalar>& cWeights,
+    List<label>& decomp
+)
+{
+    // Default setup
+    enum configs kahipConfig = configs::FAST;
+    double imbalance = 0.01;
+    int seed = 0;
+    bool verbose = false;
+
+    const dictionary* coeffsDictPtr =
+        decompositionDict_.subDictPtr("kahipCoeffs");
+
+    #if WM_LABEL_SIZE == 64
+    if (xadj.size()-1 > INT_MAX)
+    {
+        FatalErrorInFunction
+            << "Cannot decompose " << (xadj.size()-1) << " cells," << nl
+            << "Exceeded integer limit of " << INT_MAX << nl
+            << exit(FatalError);
+    }
+    #endif
+
+    int numCells = xadj.size()-1;
+
+    // Cell weights (so on the vertices of the dual)
+    List<int> cellWeights;
+
+    // Check for externally provided cellweights and if so initialise weights
+    const scalar minWeights = gMin(cWeights);
+    if (!cWeights.empty())
+    {
+        if (minWeights <= 0)
+        {
+            WarningInFunction
+                << "Illegal minimum weight " << minWeights
+                << endl;
+        }
+
+        if (cWeights.size() != numCells)
+        {
+            FatalErrorInFunction
+                << "Number of cell weights " << cWeights.size()
+                << " does not equal number of cells " << numCells
+                << exit(FatalError);
+        }
+
+        // Convert to integers.
+        cellWeights.setSize(cWeights.size());
+        forAll(cellWeights, i)
+        {
+            cellWeights[i] = int(cWeights[i]/minWeights);
+        }
+    }
+
+    // Additional sizing parameters (testing only)
+    std::map<std::string, std::vector<int>> sizingParams;
+
+    // Check for user supplied weights and decomp options
+    if (coeffsDictPtr)
+    {
+        const dictionary& coeffDict = *coeffsDictPtr;
+
+        //- Find the key in the dictionary and return the corresponding
+        //  enumeration element based on its name.
+        //  Return the default value if the key was not found in the dictionary.
+        //  Fatal if the enumerated name was incorrect.
+
+        kahipConfig =
+            configNames.lookupOrDefault("config", coeffDict, kahipConfig);
+
+        coeffDict.readIfPresent("imbalance", imbalance);
+        coeffDict.readIfPresent("verbose", verbose);
+
+        Info<< "kahipDecomp :"
+            << " config=" << configNames[kahipConfig]
+            << " imbalance=" << imbalance;
+
+        List<int> labels;
+        if
+        (
+            coeffDict.readIfPresent("hierarchy", labels)
+         && !labels.empty()
+        )
+        {
+            std::vector<int> vec;
+            vec.reserve(labels.size()+1);
+
+            // Verify sizing
+
+            int n = 1;
+            for (auto val : labels)
+            {
+                n *= val;
+                vec.push_back(val);
+            }
+
+            if (n != nProcessors_)
+            {
+                // Size mismatch. Try to correct.
+
+                if (nProcessors_ % n)
+                {
+                    WarningInFunction
+                        << "Mismatch in number of processors and "
+                        << "hierarchy specified" << flatOutput(labels) << endl;
+
+                    vec.clear();
+                }
+                else
+                {
+                    // Evenly divisible, add extra hierarchy level
+                    vec.push_back(nProcessors_ / n);
+                }
+            }
+
+            if (!vec.empty())
+            {
+                sizingParams["hierarchy"] = std::move(vec);
+
+                Info<< " hierarchy=" << flatOutput(labels);
+            }
+        }
+
+        if
+        (
+            coeffDict.readIfPresent("distance", labels)
+         && !labels.empty()
+        )
+        {
+            std::vector<int> vec(labels.size());
+
+            forAll(labels, i)
+            {
+                vec[i] = labels[i];
+            }
+
+            sizingParams["distance"] = std::move(vec);
+
+            Info<< " distance=" << flatOutput(labels);
+        }
+
+        if (coeffDict.readIfPresent("seed", seed))
+        {
+            Info<< " seed=" << seed;
+        }
+
+        Info<< endl;
+    }
+    else
+    {
+        Info<< "kahipDecomp :"
+            << " config=" << configNames[kahipConfig]
+            << " imbalance=" << imbalance << endl;
+    }
+
+    // Number of partitions
+    int nParts = nProcessors_;
+
+    // Output: number of cut edges
+    int edgeCut = 0;
+
+    #if WM_LABEL_SIZE == 32
+
+    // Input:
+    int* xadjPtr   = const_cast<UList<int>&>(xadj).begin();
+    int* adjncyPtr = const_cast<UList<int>&>(adjncy).begin();
+
+    // Output: cell -> processor addressing
+    decomp.setSize(numCells);
+    int* decompPtr = decomp.begin();
+
+    #elif WM_LABEL_SIZE == 64
+
+    // input (copy)
+    List<int> xadjCopy(xadj.size());
+    List<int> adjncyCopy(adjncy.size());
+
+    forAll(xadj,i)
+    {
+        xadjCopy[i] = xadj[i];
+    }
+    forAll(adjncy,i)
+    {
+        adjncyCopy[i] = adjncy[i];
+    }
+
+    int* xadjPtr   = xadjCopy.begin();
+    int* adjncyPtr = adjncyCopy.begin();
+
+    if (decomp.size() != numCells)
+    {
+        decomp.clear();
+    }
+
+    // Output: cell -> processor addressing
+    List<int> decompCopy(numCells);
+    int* decompPtr = decompCopy.begin();
+    #endif
+
+#if 0 // WIP: #ifdef KAFFPA_CPP_INTERFACE
+    kaffpa_cpp
+    (
+        &numCells,          // num vertices in graph
+        (cellWeights.size() ? cellWeights.begin() : nullptr), // vertex wts
+        xadjPtr,            // indexing into adjncy
+        nullptr,            // edge wts
+        adjncyPtr,          // neighbour info
+        &nParts,            // nparts
+        &imbalance,         // amount of imbalance allowed
+        !verbose,           // suppress output
+        seed,               // for random
+        int(kahipConfig),
+        &edgeCut,           // [output]
+        decompPtr,          // [output]
+        sizingParams
+    );
+#else
+    kaffpa
+    (
+        &numCells,          // num vertices in graph
+        (cellWeights.size() ? cellWeights.begin() : nullptr), // vertex wts
+        xadjPtr,            // indexing into adjncy
+        nullptr,            // edge wts
+        adjncyPtr,          // neighbour info
+        &nParts,            // nparts
+        &imbalance,         // amount of imbalance allowed
+        !verbose,           // suppress output
+        seed,               // for random
+        int(kahipConfig),
+        &edgeCut,           // [output]
+        decompPtr           // [output]
+    );
+#endif
+
+    #if WM_LABEL_SIZE == 64
+
+    // Drop input copy
+    xadjCopy.clear();
+    adjncyCopy.clear();
+
+    // Copy back to List<label>
+    decomp.setSize(numCells);
+    forAll(decompCopy, i)
+    {
+        decomp[i] = decompCopy[i];
+    }
+
+    decompCopy.clear();
+    #endif
+
+    return edgeCut;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::kahipDecomp::kahipDecomp(const dictionary& decompositionDict)
+:
+    metisLikeDecomp(decompositionDict)
+{}
+
+
+// ************************************************************************* //

src/parallel/decompose/kahipDecomp/kahipDecomp.H

@@ -0,0 +1,143 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     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::kahipDecomp
+
+Description
+    Domain decomposition using KaHIP
+    http://algo2.iti.kit.edu/documents/kahip/
+
+    When run in parallel will collect the entire graph on to the master,
+    decompose and send back.
+
+    \verbatim
+    numberOfSubdomains   N;
+    method               kahip;
+
+    kahipCoeffs
+    {
+        config          fast;
+        imbalance       0.01;
+    }
+    \endverbatim
+
+    Where the coefficients dictionary is optional, as are all its entries:
+    \table
+        Property  | Description                                 | Default value
+        config    | fast / eco / strong                         | fast
+        imbalance | imbalance on cells between domains          | 0.01
+        seed      | initial value for random number generator   | 0
+    \endtable
+
+SourceFiles
+    kahipDecomp.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef kahipDecomp_H
+#define kahipDecomp_H
+
+#include "metisLikeDecomp.H"
+#include "Enum.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                         Class kahipDecomp Declaration
+\*---------------------------------------------------------------------------*/
+
+class kahipDecomp
+:
+    public metisLikeDecomp
+{
+
+    // Private Member Functions
+
+        //- Call kahip with options from dictionary.
+        virtual label decomposeSerial
+        (
+            const UList<label>& adjncy,
+            const UList<label>& xadj,
+            const UList<scalar>& cellWeights,
+            List<label>& decomp
+        );
+
+        //- Disallow default bitwise copy construct and assignment
+        void operator=(const kahipDecomp&) = delete;
+        kahipDecomp(const kahipDecomp&) = delete;
+
+
+public:
+
+    //- The predefined KaHIP configuration types
+    enum class configs
+    {
+        FAST = 0, //!< default
+        ECO = 1,
+        STRONG = 2,
+        FASTSOCIAL = 3,
+        ECOSOCIAL = 4,
+        STRONGSOCIAL = 5,
+    };
+
+
+    //- The selection names for predefined KaHIP configurations
+    static const Enum<configs> configNames;
+
+
+    //- Runtime type information
+    TypeName("kahip");
+
+
+    // Constructors
+
+        //- Construct given the decomposition dictionary
+        kahipDecomp(const dictionary& decompositionDict);
+
+
+    //- Destructor
+    virtual ~kahipDecomp()
+    {}
+
+
+    // Member Functions
+
+        virtual bool parallelAware() const
+        {
+            return true;
+        }
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/parallel/decompose/metisDecomp/metisDecomp.C

@@ -45,13 +45,12 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-Foam::label Foam::metisDecomp::decompose
+Foam::label Foam::metisDecomp::decomposeSerial
 (
-    const List<label>& adjncy,
+    const UList<label>& adjncy,
-    const List<label>& xadj,
+    const UList<label>& xadj,
-    const scalarField& cWeights,
+    const UList<scalar>& cWeights,
-
+    List<label>& decomp
-    List<label>& finalDecomp
 )
 {
     // Method of decomposition
@@ -59,6 +58,9 @@ Foam::label Foam::metisDecomp::decompose
     // k-way: multi-level k-way
     word method("recursive");
 
+    const dictionary* coeffsDictPtr =
+        decompositionDict_.subDictPtr("metisCoeffs");
+
     label numCells = xadj.size()-1;
 
     // Decomposition options
@@ -77,7 +79,7 @@ Foam::label Foam::metisDecomp::decompose
 
 
     // Check for externally provided cellweights and if so initialise weights
-    scalar minWeights = gMin(cWeights);
+    const scalar minWeights = gMin(cWeights);
     if (cWeights.size() > 0)
     {
         if (minWeights <= 0)
@@ -105,14 +107,13 @@ Foam::label Foam::metisDecomp::decompose
 
 
     // Check for user supplied weights and decomp options
-    if (decompositionDict_.found("metisCoeffs"))
+    if (coeffsDictPtr)
     {
-        const dictionary& metisCoeffs =
+        const dictionary& coeffDict = *coeffsDictPtr;
-            decompositionDict_.subDict("metisCoeffs");
 
         word weightsFile;
 
-        if (metisCoeffs.readIfPresent("method", method))
+        if (coeffDict.readIfPresent("method", method))
         {
             if (method != "recursive" && method != "k-way")
             {
@@ -127,7 +128,7 @@ Foam::label Foam::metisDecomp::decompose
                 << nl << endl;
         }
 
-        if (metisCoeffs.readIfPresent("options", options))
+        if (coeffDict.readIfPresent("options", options))
         {
             if (options.size() != METIS_NOPTIONS)
             {
@@ -142,7 +143,7 @@ Foam::label Foam::metisDecomp::decompose
                 << nl << endl;
         }
 
-        if (metisCoeffs.readIfPresent("processorWeights", processorWeights))
+        if (coeffDict.readIfPresent("processorWeights", processorWeights))
         {
             processorWeights /= sum(processorWeights);
 
@@ -161,7 +162,7 @@ Foam::label Foam::metisDecomp::decompose
     label nProcs = nProcessors_;
 
     // Output: cell -> processor addressing
-    finalDecomp.setSize(numCells);
+    decomp.setSize(numCells);
 
     // Output: number of cut edges
     label edgeCut = 0;
@@ -172,17 +173,17 @@ Foam::label Foam::metisDecomp::decompose
         (
             &numCells,          // num vertices in graph
             &ncon,              // num balancing constraints
-            const_cast<List<label>&>(xadj).begin(),   // indexing into adjncy
+            const_cast<UList<label>&>(xadj).begin(),   // indexing into adjncy
-            const_cast<List<label>&>(adjncy).begin(), // neighbour info
+            const_cast<UList<label>&>(adjncy).begin(), // neighbour info
-            cellWeights.begin(),// vertexweights
+            cellWeights.begin(),// vertex wts
             nullptr,               // vsize: total communication vol
-            faceWeights.begin(),// edgeweights
+            faceWeights.begin(),// edge wts
             &nProcs,            // nParts
             processorWeights.begin(),   // tpwgts
             nullptr,               // ubvec: processor imbalance (default)
             options.begin(),
             &edgeCut,
-            finalDecomp.begin()
+            decomp.begin()
         );
     }
     else
@@ -191,17 +192,17 @@ Foam::label Foam::metisDecomp::decompose
         (
             &numCells,          // num vertices in graph
             &ncon,              // num balancing constraints
-            const_cast<List<label>&>(xadj).begin(),   // indexing into adjncy
+            const_cast<UList<label>&>(xadj).begin(),   // indexing into adjncy
-            const_cast<List<label>&>(adjncy).begin(), // neighbour info
+            const_cast<UList<label>&>(adjncy).begin(), // neighbour info
-            cellWeights.begin(),// vertexweights
+            cellWeights.begin(),// vertex wts
             nullptr,               // vsize: total communication vol
-            faceWeights.begin(),// edgeweights
+            faceWeights.begin(),// edge wts
             &nProcs,            // nParts
             processorWeights.begin(),   // tpwgts
             nullptr,               // ubvec: processor imbalance (default)
             options.begin(),
             &edgeCut,
-            finalDecomp.begin()
+            decomp.begin()
         );
     }
 
@@ -213,116 +214,8 @@ Foam::label Foam::metisDecomp::decompose
 
 Foam::metisDecomp::metisDecomp(const dictionary& decompositionDict)
 :
-    decompositionMethod(decompositionDict)
+    metisLikeDecomp(decompositionDict)
 {}
 
 
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::labelList Foam::metisDecomp::decompose
-(
-    const polyMesh& mesh,
-    const pointField& points,
-    const scalarField& pointWeights
-)
-{
-    if (points.size() != mesh.nCells())
-    {
-        FatalErrorInFunction
-            << "Can use this decomposition method only for the whole mesh"
-            << endl
-            << "and supply one coordinate (cellCentre) for every cell." << endl
-            << "The number of coordinates " << points.size() << endl
-            << "The number of cells in the mesh " << mesh.nCells()
-            << exit(FatalError);
-    }
-
-    CompactListList<label> cellCells;
-    calcCellCells
-    (
-        mesh,
-        identity(mesh.nCells()),
-        mesh.nCells(),
-        false,
-        cellCells
-    );
-
-    // Decompose using default weights
-    labelList decomp;
-    decompose(cellCells.m(), cellCells.offsets(), pointWeights, decomp);
-
-    return decomp;
-}
-
-
-Foam::labelList Foam::metisDecomp::decompose
-(
-    const polyMesh& mesh,
-    const labelList& agglom,
-    const pointField& agglomPoints,
-    const scalarField& agglomWeights
-)
-{
-    if (agglom.size() != mesh.nCells())
-    {
-        FatalErrorInFunction
-            << "Size of cell-to-coarse map " << agglom.size()
-            << " differs from number of cells in mesh " << mesh.nCells()
-            << exit(FatalError);
-    }
-
-    // Make Metis CSR (Compressed Storage Format) storage
-    //   adjncy      : contains neighbours (= edges in graph)
-    //   xadj(celli) : start of information in adjncy for celli
-
-    CompactListList<label> cellCells;
-    calcCellCells(mesh, agglom, agglomPoints.size(), false, cellCells);
-
-    // Decompose using default weights
-    labelList finalDecomp;
-    decompose(cellCells.m(), cellCells.offsets(), agglomWeights, finalDecomp);
-
-
-    // Rework back into decomposition for original mesh
-    labelList fineDistribution(agglom.size());
-
-    forAll(fineDistribution, i)
-    {
-        fineDistribution[i] = finalDecomp[agglom[i]];
-    }
-
-    return fineDistribution;
-}
-
-
-Foam::labelList Foam::metisDecomp::decompose
-(
-    const labelListList& globalCellCells,
-    const pointField& cellCentres,
-    const scalarField& cellWeights
-)
-{
-    if (cellCentres.size() != globalCellCells.size())
-    {
-        FatalErrorInFunction
-            << "Inconsistent number of cells (" << globalCellCells.size()
-            << ") and number of cell centres (" << cellCentres.size()
-            << ")." << exit(FatalError);
-    }
-
-
-    // Make Metis CSR (Compressed Storage Format) storage
-    //   adjncy      : contains neighbours (= edges in graph)
-    //   xadj(celli) : start of information in adjncy for celli
-
-    CompactListList<label> cellCells(globalCellCells);
-
-    // Decompose using default weights
-    labelList decomp;
-    decompose(cellCells.m(), cellCells.offsets(), cellWeights, decomp);
-
-    return decomp;
-}
-
-
 // ************************************************************************* //

src/parallel/decompose/metisDecomp/metisDecomp.H

@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
-     \\/     M anipulation  |
+     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -27,6 +27,29 @@ Class
 Description
     Metis domain decomposition
 
+    When run in parallel will collect the entire graph on to the master,
+    decompose and send back.
+
+    \verbatim
+    numberOfSubdomains   N;
+    method               metis;
+
+    metisCoeffs
+    {
+        method          recursive;  // k-way
+        options          ( ...);
+        processorWeights (   ... );
+    }
+    \endverbatim
+
+    Where the coefficients dictionary is optional, as are all its entries:
+    \table
+        Property  | Description                             | Default value
+        method    | recursive / k-way                       | recursive
+        options   | metis options                           |
+        processorWeights | list of weighting per partition  |
+    \endtable
+
 SourceFiles
     metisDecomp.C
 
@@ -35,7 +58,7 @@ SourceFiles
 #ifndef metisDecomp_H
 #define metisDecomp_H
 
-#include "decompositionMethod.H"
+#include "metisLikeDecomp.H"
 
 namespace Foam
 {
@@ -46,23 +69,23 @@ namespace Foam
 
 class metisDecomp
 :
-    public decompositionMethod
+    public metisLikeDecomp
 {
 
     // Private Member Functions
 
         //- Call Metis with options from dictionary.
-        label decompose
+        virtual label decomposeSerial
         (
-            const List<label>& adjncy,
+            const UList<label>& adjncy,
-            const List<label>& xadj,
+            const UList<label>& xadj,
-            const scalarField& cellWeights,
+            const UList<scalar>& cellWeights,
-            List<label>& finalDecomp
+            List<label>& decomp
         );
 
         //- Disallow default bitwise copy construct and assignment
-        void operator=(const metisDecomp&);
+        void operator=(const metisDecomp&) = delete;
-        metisDecomp(const metisDecomp&);
+        metisDecomp(const metisDecomp&) = delete;
 
 
 public:
@@ -74,7 +97,7 @@ public:
     // Constructors
 
         //- Construct given the decomposition dictionary
-        metisDecomp(const dictionary&);
+        metisDecomp(const dictionary& decompositionDict);
 
 
     //- Destructor
@@ -86,52 +109,9 @@ public:
 
         virtual bool parallelAware() const
         {
-            // Metis does not know about proc boundaries
+            return true;
-            return false;
         }
 
-        //- Inherit decompose from decompositionMethod
-        using decompositionMethod::decompose;
-
-        //- Return for every coordinate the wanted processor number. Use the
-        //  mesh connectivity (if needed)
-        //  Weights get normalised so the minimum value is 1 before truncation
-        //  to an integer so the weights should be multiples of the minimum
-        //  value. The overall sum of weights might otherwise overflow.
-        virtual labelList decompose
-        (
-            const polyMesh& mesh,
-            const pointField& points,
-            const scalarField& pointWeights
-        );
-
-        //- Return for every coordinate the wanted processor number. Gets
-        //  passed agglomeration map (from fine to coarse cells) and coarse cell
-        //  location. Can be overridden by decomposers that provide this
-        //  functionality natively.
-        //  See note on weights above.
-        virtual labelList decompose
-        (
-            const polyMesh& mesh,
-            const labelList& agglom,
-            const pointField& regionPoints,
-            const scalarField& regionWeights
-        );
-
-        //- Return for every coordinate the wanted processor number. Explicitly
-        //  provided mesh connectivity.
-        //  The connectivity is equal to mesh.cellCells() except for
-        //  - in parallel the cell numbers are global cell numbers (starting
-        //    from 0 at processor0 and then incrementing all through the
-        //    processors)
-        //  - the connections are across coupled patches
-        //  See note on weights above.
-        virtual labelList decompose
-        (
-            const labelListList& globalCellCells,
-            const pointField& cc,
-            const scalarField& cWeights
-        );
 };
 
 

src/parallel/decompose/ptscotchDecomp/ptscotchDecomp.C

@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2015-2017 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -228,12 +228,17 @@ License
 namespace Foam
 {
     defineTypeNameAndDebug(ptscotchDecomp, 0);
-
     addToRunTimeSelectionTable(decompositionMethod, ptscotchDecomp, dictionary);
 }
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+void Foam::ptscotchDecomp::graphPath(const polyMesh& mesh)
+{
+    graphPath_ = mesh.time().path()/mesh.name();
+}
+
+
 void Foam::ptscotchDecomp::check(const int retVal, const char* str)
 {
     if (retVal)
@@ -248,11 +253,9 @@ void Foam::ptscotchDecomp::check(const int retVal, const char* str)
 ////- Does prevention of 0 cell domains and calls ptscotch.
 //Foam::label Foam::ptscotchDecomp::decomposeZeroDomains
 //(
-//    const fileName& meshPath,
+//    const UList<label>& initadjncy,
-//    const List<label>& initadjncy,
+//    const UList<label>& initxadj,
-//    const List<label>& initxadj,
+//    const UList<scalar>& initcWeights,
-//    const scalarField& initcWeights,
-//
 //    List<label>& finalDecomp
 //) const
 //{
@@ -272,7 +275,6 @@ void Foam::ptscotchDecomp::check(const int retVal, const char* str)
 //    {
 //        return decompose
 //        (
-//            meshPath,
 //            initadjncy,
 //            initxadj,
 //            initcWeights,
@@ -381,7 +383,7 @@ void Foam::ptscotchDecomp::check(const int retVal, const char* str)
 //
 //
 //    // Do decomposition as normal. Sets finalDecomp.
-//    label result = decompose(meshPath, adjncy, xadj, cWeights, finalDecomp);
+//    label result = decompose(adjncy, xadj, cWeights, finalDecomp);
 //
 //
 //    if (debug)
@@ -437,23 +439,19 @@ void Foam::ptscotchDecomp::check(const int retVal, const char* str)
 //}
 
 
-// Call scotch with options from dictionary.
 Foam::label Foam::ptscotchDecomp::decompose
 (
-    const fileName& meshPath,
+    const UList<label>& adjncy,
-    const List<label>& adjncy,
+    const UList<label>& xadj,
-    const List<label>& xadj,
+    const UList<scalar>& cWeights,
-    const scalarField& cWeights,
     List<label>& finalDecomp
 ) const
 {
-    List<label> dummyAdjncy(1);
+    List<label> dummyAdjncy { 0 };
-    List<label> dummyXadj(1);
+    List<label> dummyXadj { 0 };
-    dummyXadj[0] = 0;
 
     return decompose
     (
-        meshPath,
         adjncy.size(),
         (adjncy.size() ? adjncy.begin() : dummyAdjncy.begin()),
         xadj.size(),
@@ -464,16 +462,13 @@ Foam::label Foam::ptscotchDecomp::decompose
 }
 
 
-// Call scotch with options from dictionary.
 Foam::label Foam::ptscotchDecomp::decompose
 (
-    const fileName& meshPath,
     const label adjncySize,
     const label adjncy[],
     const label xadjSize,
     const label xadj[],
-    const scalarField& cWeights,
+    const UList<scalar>& cWeights,
-
     List<label>& finalDecomp
 ) const
 {
@@ -492,7 +487,7 @@ Foam::label Foam::ptscotchDecomp::decompose
         {
             OFstream str
             (
-               meshPath + "_" + Foam::name(Pstream::myProcNo()) + ".dgr"
+               graphPath_ + "_" + Foam::name(Pstream::myProcNo()) + ".dgr"
             );
 
             Pout<< "Dumping Scotch graph file to " << str.name() << endl
@@ -572,8 +567,8 @@ Foam::label Foam::ptscotchDecomp::decompose
 
     // Check for externally provided cellweights and if so initialise weights
 
-    scalar minWeights = gMin(cWeights);
+    const scalar minWeights = gMin(cWeights);
-    scalar maxWeights = gMax(cWeights);
+    const scalar maxWeights = gMax(cWeights);
 
     if (maxWeights > minWeights)
     {
@@ -719,7 +714,10 @@ Foam::label Foam::ptscotchDecomp::decompose
         }
         check
         (
-            SCOTCH_archCmpltw(&archdat, nProcessors_, processorWeights.begin()),
+            SCOTCH_archCmpltw
+            (
+                &archdat, nProcessors_, processorWeights.begin()
+            ),
             "SCOTCH_archCmpltw"
         );
     }
@@ -819,6 +817,9 @@ Foam::labelList Foam::ptscotchDecomp::decompose
     const scalarField& pointWeights
 )
 {
+    // Where to write graph
+    graphPath(mesh);
+
     if (points.size() != mesh.nCells())
     {
         FatalErrorInFunction
@@ -835,7 +836,6 @@ Foam::labelList Foam::ptscotchDecomp::decompose
     //   adjncy      : contains neighbours (= edges in graph)
     //   xadj(celli) : start of information in adjncy for celli
 
-
     CompactListList<label> cellCells;
     calcCellCells
     (
@@ -847,22 +847,15 @@ Foam::labelList Foam::ptscotchDecomp::decompose
     );
 
     // Decompose using default weights
-    List<label> finalDecomp;
+    labelList decomp;
     decompose
     (
-        mesh.time().path()/mesh.name(),
         cellCells.m(),
         cellCells.offsets(),
         pointWeights,
-        finalDecomp
+        decomp
     );
 
-    // Copy back to labelList
-    labelList decomp(points.size());
-    forAll(decomp, i)
-    {
-        decomp[i] = finalDecomp[i];
-    }
     return decomp;
 }
 
@@ -875,6 +868,9 @@ Foam::labelList Foam::ptscotchDecomp::decompose
     const scalarField& pointWeights
 )
 {
+    // Where to write graph
+    graphPath(mesh);
+
     if (agglom.size() != mesh.nCells())
     {
         FatalErrorInFunction
@@ -898,14 +894,13 @@ Foam::labelList Foam::ptscotchDecomp::decompose
     );
 
     // Decompose using weights
-    List<label> finalDecomp;
+    labelList decomp;
     decompose
     (
-        mesh.time().path()/mesh.name(),
         cellCells.m(),
         cellCells.offsets(),
         pointWeights,
-        finalDecomp
+        decomp
     );
 
     // Rework back into decomposition for original mesh
@@ -913,7 +908,7 @@ Foam::labelList Foam::ptscotchDecomp::decompose
 
     forAll(fineDistribution, i)
     {
-        fineDistribution[i] = finalDecomp[agglom[i]];
+        fineDistribution[i] = decomp[agglom[i]];
     }
 
     return fineDistribution;
@@ -927,6 +922,9 @@ Foam::labelList Foam::ptscotchDecomp::decompose
     const scalarField& cWeights
 )
 {
+    // Where to write graph
+    graphPath_ = "ptscotch";
+
     if (cellCentres.size() != globalCellCells.size())
     {
         FatalErrorInFunction
@@ -943,22 +941,15 @@ Foam::labelList Foam::ptscotchDecomp::decompose
     CompactListList<label> cellCells(globalCellCells);
 
     // Decompose using weights
-    List<label> finalDecomp;
+    labelList decomp;
     decompose
     (
-        "ptscotch",
         cellCells.m(),
         cellCells.offsets(),
         cWeights,
-        finalDecomp
+        decomp
     );
 
-    // Copy back to labelList
-    labelList decomp(cellCentres.size());
-    forAll(decomp, i)
-    {
-        decomp[i] = finalDecomp[i];
-    }
     return decomp;
 }
 

src/parallel/decompose/ptscotchDecomp/ptscotchDecomp.H

@@ -64,36 +64,43 @@ class ptscotchDecomp
 :
     public decompositionMethod
 {
+    // Private data
+
+        //- Output path and name for optional grf file.
+        fileName graphPath_;
+
+
     // Private Member Functions
 
+        //- Set graph path and name
+        void graphPath(const polyMesh& mesh);
+
         //- Check and print error message
         static void check(const int, const char*);
 
         //- Decompose. Handles size 0 arrays
         label decompose
         (
-            const fileName& meshPath,
+            const UList<label>& adjncy,
-            const List<label>& adjncy,
+            const UList<label>& xadj,
-            const List<label>& xadj,
+            const UList<scalar>& cWeights,
-            const scalarField& cWeights,
             List<label>& finalDecomp
         ) const;
 
         //- Low level decompose
         label decompose
         (
-            const fileName& meshPath,
             const label adjncySize,
             const label adjncy[],
             const label xadjSize,
             const label xadj[],
-            const scalarField& cWeights,
+            const UList<scalar>& cWeights,
             List<label>& finalDecomp
         ) const;
 
         //- Disallow default bitwise copy construct and assignment
-        void operator=(const ptscotchDecomp&);
+        void operator=(const ptscotchDecomp&) = delete;
-        ptscotchDecomp(const ptscotchDecomp&);
+        ptscotchDecomp(const ptscotchDecomp&) = delete;
 
 
 public:

src/parallel/decompose/scotchDecomp/scotchDecomp.C

@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2015-2017 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -126,8 +126,6 @@ License
 #include "floatScalar.H"
 #include "Time.H"
 #include "OFstream.H"
-#include "globalIndex.H"
-#include "SubField.H"
 
 extern "C"
 {
@@ -161,6 +159,12 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
+void Foam::scotchDecomp::graphPath(const polyMesh& mesh)
+{
+    graphPath_ = mesh.time().path()/mesh.name() + ".grf";
+}
+
+
 void Foam::scotchDecomp::check(const int retVal, const char* str)
 {
     if (retVal)
@@ -172,185 +176,52 @@ void Foam::scotchDecomp::check(const int retVal, const char* str)
 }
 
 
-Foam::label Foam::scotchDecomp::decompose
+Foam::label Foam::scotchDecomp::decomposeSerial
 (
-    const fileName& meshPath,
+    const UList<label>& adjncy,
-    const List<label>& adjncy,
+    const UList<label>& xadj,
-    const List<label>& xadj,
+    const UList<scalar>& cWeights,
-    const scalarField& cWeights,
+    List<label>& decomp
-
-    List<label>& finalDecomp
 )
 {
-    if (!Pstream::parRun())
+    const dictionary* coeffsDictPtr =
-    {
+        decompositionDict_.subDictPtr("scotchCoeffs");
-        decomposeOneProc
-        (
-            meshPath,
-            adjncy,
-            xadj,
-            cWeights,
-            finalDecomp
-        );
-    }
-    else
-    {
-        if (debug)
-        {
-            Info<< "scotchDecomp : running in parallel."
-                << " Decomposing all of graph on master processor." << endl;
-        }
-        globalIndex globalCells(xadj.size()-1);
-        label nTotalConnections = returnReduce(adjncy.size(), sumOp<label>());
 
-        // Send all to master. Use scheduled to save some storage.
-        if (Pstream::master())
-        {
-            Field<label> allAdjncy(nTotalConnections);
-            Field<label> allXadj(globalCells.size()+1);
-            scalarField allWeights(globalCells.size());
-
-            // Insert my own
-            label nTotalCells = 0;
-            forAll(cWeights, celli)
-            {
-                allXadj[nTotalCells] = xadj[celli];
-                allWeights[nTotalCells++] = cWeights[celli];
-            }
-            nTotalConnections = 0;
-            forAll(adjncy, i)
-            {
-                allAdjncy[nTotalConnections++] = adjncy[i];
-            }
-
-            for (int slave=1; slave<Pstream::nProcs(); slave++)
-            {
-                IPstream fromSlave(Pstream::commsTypes::scheduled, slave);
-                Field<label> nbrAdjncy(fromSlave);
-                Field<label> nbrXadj(fromSlave);
-                scalarField nbrWeights(fromSlave);
-
-                // Append.
-                //label procStart = nTotalCells;
-                forAll(nbrXadj, celli)
-                {
-                    allXadj[nTotalCells] = nTotalConnections+nbrXadj[celli];
-                    allWeights[nTotalCells++] = nbrWeights[celli];
-                }
-                // No need to renumber xadj since already global.
-                forAll(nbrAdjncy, i)
-                {
-                    allAdjncy[nTotalConnections++] = nbrAdjncy[i];
-                }
-            }
-            allXadj[nTotalCells] = nTotalConnections;
-
-
-            Field<label> allFinalDecomp;
-            decomposeOneProc
-            (
-                meshPath,
-                allAdjncy,
-                allXadj,
-                allWeights,
-                allFinalDecomp
-            );
-
-
-            // Send allFinalDecomp back
-            for (int slave=1; slave<Pstream::nProcs(); slave++)
-            {
-                OPstream toSlave(Pstream::commsTypes::scheduled, slave);
-                toSlave << SubField<label>
-                (
-                    allFinalDecomp,
-                    globalCells.localSize(slave),
-                    globalCells.offset(slave)
-                );
-            }
-            // Get my own part (always first)
-            finalDecomp = SubField<label>
-            (
-                allFinalDecomp,
-                globalCells.localSize()
-            );
-        }
-        else
-        {
-            // Send my part of the graph (already in global numbering)
-            {
-                OPstream toMaster
-                (
-                    Pstream::commsTypes::scheduled,
-                    Pstream::masterNo()
-                );
-                toMaster<< adjncy << SubField<label>(xadj, xadj.size()-1)
-                    << cWeights;
-            }
-
-            // Receive back decomposition
-            IPstream fromMaster
-            (
-                Pstream::commsTypes::scheduled,
-                Pstream::masterNo()
-            );
-            fromMaster >> finalDecomp;
-        }
-    }
-    return 0;
-}
-
-
-// Call scotch with options from dictionary.
-Foam::label Foam::scotchDecomp::decomposeOneProc
-(
-    const fileName& meshPath,
-    const List<label>& adjncy,
-    const List<label>& xadj,
-    const scalarField& cWeights,
-
-    List<label>& finalDecomp
-)
-{
     // Dump graph
-    if (decompositionDict_.found("scotchCoeffs"))
+    if (coeffsDictPtr && coeffsDictPtr->lookupOrDefault("writeGraph", false))
     {
-        const dictionary& scotchCoeffs =
+        OFstream str(graphPath_);
-            decompositionDict_.subDict("scotchCoeffs");
-
-        if (scotchCoeffs.lookupOrDefault("writeGraph", false))
-        {
-            OFstream str(meshPath + ".grf");
 
         Info<< "Dumping Scotch graph file to " << str.name() << endl
             << "Use this in combination with gpart." << endl;
 
-            label version = 0;
+        const label version = 0;
         str << version << nl;
         // Numer of vertices
         str << xadj.size()-1 << ' ' << adjncy.size() << nl;
-            // Numbering starts from 0
-            label baseval = 0;
-            // Has weights?
-            label hasEdgeWeights = 0;
-            label hasVertexWeights = 0;
-            label numericflag = 10*hasEdgeWeights+hasVertexWeights;
-            str << baseval << ' ' << numericflag << nl;
-            for (label celli = 0; celli < xadj.size()-1; celli++)
-            {
-                label start = xadj[celli];
-                label end = xadj[celli+1];
-                str << end-start;
 
-                for (label i = start; i < end; i++)
+        // Numbering starts from 0
+        const label baseval = 0;
+        // Has weights?
+        const label hasEdgeWeights = 0;
+        const label hasVertexWeights = 0;
+        const label numericflag = 10*hasEdgeWeights+hasVertexWeights;
+        str << baseval << ' ' << numericflag << nl;
+
+        for (label celli = 0; celli < xadj.size()-1; ++celli)
+        {
+            const label start = xadj[celli];
+            const label end = xadj[celli+1];
+
+            str << end-start; // size
+
+            for (label i = start; i < end; ++i)
             {
                 str << ' ' << adjncy[i];
             }
             str << nl;
         }
     }
-    }
-
 
     // Strategy
     // ~~~~~~~~
@@ -359,13 +230,10 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
     SCOTCH_Strat stradat;
     check(SCOTCH_stratInit(&stradat), "SCOTCH_stratInit");
 
-    if (decompositionDict_.found("scotchCoeffs"))
+    if (coeffsDictPtr)
     {
-        const dictionary& scotchCoeffs =
-            decompositionDict_.subDict("scotchCoeffs");
-
         string strategy;
-        if (scotchCoeffs.readIfPresent("strategy", strategy))
+        if (coeffsDictPtr->readIfPresent("strategy", strategy))
         {
             if (debug)
             {
@@ -384,10 +252,9 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
 
     List<label> velotab;
 
-
     // Check for externally provided cellweights and if so initialise weights
     // Note: min, not gMin since routine runs on master only.
-    scalar minWeights = min(cWeights);
+    const scalar minWeights = min(cWeights);
     if (!cWeights.empty())
     {
         if (minWeights <= 0)
@@ -431,7 +298,6 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
     }
 
 
-
     SCOTCH_Graph grafdat;
     check(SCOTCH_graphInit(&grafdat), "SCOTCH_graphInit");
     check
@@ -462,14 +328,12 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
     check(SCOTCH_archInit(&archdat), "SCOTCH_archInit");
 
     List<label> processorWeights;
-    if (decompositionDict_.found("scotchCoeffs"))
+    if
-    {
+    (
-        const dictionary& scotchCoeffs =
+        coeffsDictPtr
-            decompositionDict_.subDict("scotchCoeffs");
+     && coeffsDictPtr->readIfPresent("processorWeights", processorWeights)
-
+     && processorWeights.size()
-        scotchCoeffs.readIfPresent("processorWeights", processorWeights);
+    )
-    }
-    if (processorWeights.size())
     {
         if (debug)
         {
@@ -478,7 +342,10 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
         }
         check
         (
-            SCOTCH_archCmpltw(&archdat, nProcessors_, processorWeights.begin()),
+            SCOTCH_archCmpltw
+            (
+                &archdat, nProcessors_, processorWeights.begin()
+            ),
             "SCOTCH_archCmpltw"
         );
     }
@@ -491,7 +358,7 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
         );
 
 
-        //- Hack to test clustering. Note that finalDecomp is non-compact
+        //- Hack to test clustering. Note that decomp is non-compact
         //  numbers!
         //
         ////- Set up variable sizes architecture
@@ -542,8 +409,8 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
     );
     #endif
 
-    finalDecomp.setSize(xadj.size()-1);
+    decomp.setSize(xadj.size()-1);
-    finalDecomp = 0;
+    decomp = 0;
     check
     (
         SCOTCH_graphMap
@@ -551,7 +418,7 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
             &grafdat,
             &archdat,
             &stradat,       // const SCOTCH_Strat *
-            finalDecomp.begin() // parttab
+            decomp.begin()  // parttab
         ),
         "SCOTCH_graphMap"
     );
@@ -560,9 +427,7 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
     feenableexcept(oldExcepts);
     #endif
 
-
+    //decomp.setSize(xadj.size()-1);
-
-    //finalDecomp.setSize(xadj.size()-1);
     //check
     //(
     //    SCOTCH_graphPart
@@ -570,7 +435,7 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
     //        &grafdat,
     //        nProcessors_,   // partnbr
     //        &stradat,       // const SCOTCH_Strat *
-    //        finalDecomp.begin() // parttab
+    //        decomp.begin()  // parttab
     //    ),
     //    "SCOTCH_graphPart"
     //);
@@ -590,7 +455,7 @@ Foam::label Foam::scotchDecomp::decomposeOneProc
 
 Foam::scotchDecomp::scotchDecomp(const dictionary& decompositionDict)
 :
-    decompositionMethod(decompositionDict)
+    metisLikeDecomp(decompositionDict)
 {}
 
 
@@ -603,46 +468,15 @@ Foam::labelList Foam::scotchDecomp::decompose
     const scalarField& pointWeights
 )
 {
-    if (points.size() != mesh.nCells())
+    // Where to write graph
-    {
+    graphPath(mesh);
-        FatalErrorInFunction
-            << "Can use this decomposition method only for the whole mesh"
-            << endl
-            << "and supply one coordinate (cellCentre) for every cell." << endl
-            << "The number of coordinates " << points.size() << endl
-            << "The number of cells in the mesh " << mesh.nCells()
-            << exit(FatalError);
-    }
 
-    // Calculate local or global (if Pstream::parRun()) connectivity
+    return metisLikeDecomp::decompose
-    CompactListList<label> cellCells;
-    calcCellCells
     (
         mesh,
-        identity(mesh.nCells()),
+        points,
-        mesh.nCells(),
+        pointWeights
-        true,
-        cellCells
     );
-
-    // Decompose using default weights
-    List<label> finalDecomp;
-    decompose
-    (
-        mesh.time().path()/mesh.name(),
-        cellCells.m(),
-        cellCells.offsets(),
-        pointWeights,
-        finalDecomp
-    );
-
-    // Copy back to labelList
-    labelList decomp(finalDecomp.size());
-    forAll(decomp, i)
-    {
-        decomp[i] = finalDecomp[i];
-    }
-    return decomp;
 }
 
 
@@ -654,45 +488,16 @@ Foam::labelList Foam::scotchDecomp::decompose
     const scalarField& pointWeights
 )
 {
-    if (agglom.size() != mesh.nCells())
+    // Where to write graph
-    {
+    graphPath(mesh);
-        FatalErrorInFunction
-            << "Size of cell-to-coarse map " << agglom.size()
-            << " differs from number of cells in mesh " << mesh.nCells()
-            << exit(FatalError);
-    }
 
-    // Calculate local or global (if Pstream::parRun()) connectivity
+    return metisLikeDecomp::decompose
-    CompactListList<label> cellCells;
-    calcCellCells
     (
         mesh,
         agglom,
-        agglomPoints.size(),
+        agglomPoints,
-        true,
+        pointWeights
-        cellCells
     );
-
-    // Decompose using weights
-    List<label> finalDecomp;
-    decompose
-    (
-        mesh.time().path()/mesh.name(),
-        cellCells.m(),
-        cellCells.offsets(),
-        pointWeights,
-        finalDecomp
-    );
-
-    // Rework back into decomposition for original mesh_
-    labelList fineDistribution(agglom.size());
-
-    forAll(fineDistribution, i)
-    {
-        fineDistribution[i] = finalDecomp[agglom[i]];
-    }
-
-    return fineDistribution;
 }
 
 
@@ -703,39 +508,15 @@ Foam::labelList Foam::scotchDecomp::decompose
     const scalarField& cWeights
 )
 {
-    if (cellCentres.size() != globalCellCells.size())
+    // Where to write graph
-    {
+    graphPath_ = "scotch.grf";
-        FatalErrorInFunction
-            << "Inconsistent number of cells (" << globalCellCells.size()
-            << ") and number of cell centres (" << cellCentres.size()
-            << ")." << exit(FatalError);
-    }
 
-
+    return metisLikeDecomp::decompose
-    // Make Metis CSR (Compressed Storage Format) storage
-    //   adjncy      : contains neighbours (= edges in graph)
-    //   xadj(celli) : start of information in adjncy for celli
-
-    CompactListList<label> cellCells(globalCellCells);
-
-    // Decompose using weights
-    List<label> finalDecomp;
-    decompose
     (
-        "scotch",
+        globalCellCells,
-        cellCells.m(),
+        cellCentres,
-        cellCells.offsets(),
+        cWeights
-        cWeights,
-        finalDecomp
     );
-
-    // Copy back to labelList
-    labelList decomp(finalDecomp.size());
-    forAll(decomp, i)
-    {
-        decomp[i] = finalDecomp[i];
-    }
-    return decomp;
 }
 
 

src/parallel/decompose/scotchDecomp/scotchDecomp.H

@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
-     \\/     M anipulation  |
+     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -213,7 +213,7 @@ SourceFiles
 #ifndef scotchDecomp_H
 #define scotchDecomp_H
 
-#include "decompositionMethod.H"
+#include "metisLikeDecomp.H"
 
 namespace Foam
 {
@@ -224,35 +224,34 @@ namespace Foam
 
 class scotchDecomp
 :
-    public decompositionMethod
+    public metisLikeDecomp
 {
+    // Private data
+
+        //- Output path and name for optional grf file.
+        fileName graphPath_;
+
+
     // Private Member Functions
 
+        //- Set graph path and name
+        void graphPath(const polyMesh& mesh);
+
         //- Check and print error message
         static void check(const int, const char*);
 
-        label decompose
-        (
-            const fileName& meshPath,
-            const List<label>& adjncy,
-            const List<label>& xadj,
-            const scalarField& cWeights,
-            List<label>& finalDecomp
-        );
-
         //- Decompose non-parallel
-        label decomposeOneProc
+        virtual label decomposeSerial
         (
-            const fileName& meshPath,
+            const UList<label>& adjncy,
-            const List<label>& adjncy,
+            const UList<label>& xadj,
-            const List<label>& xadj,
+            const UList<scalar>& cWeights,
-            const scalarField& cWeights,
+            List<label>& decomp
-            List<label>& finalDecomp
         );
 
         //- Disallow default bitwise copy construct and assignment
-        void operator=(const scotchDecomp&);
+        void operator=(const scotchDecomp&) = delete;
-        scotchDecomp(const scotchDecomp&);
+        scotchDecomp(const scotchDecomp&) = delete;
 
 
 public:
@@ -280,14 +279,10 @@ public:
             return true;
         }
 
-        //- Inherit decompose from decompositionMethod
+        //- Inherit all decompose methods
         using decompositionMethod::decompose;
 
-        //- Return for every coordinate the wanted processor number. Use the
+        //- Return for every coordinate the wanted processor number.
-        //  mesh connectivity (if needed)
-        //  Weights get normalised with minimum weight and truncated to
-        //  convert into integer so e.g. 3.5 is seen as 3. The overall sum
-        //  of weights might otherwise overflow.
         virtual labelList decompose
         (
             const polyMesh& mesh,
@@ -295,11 +290,7 @@ public:
             const scalarField& pointWeights
         );
 
-        //- Return for every coordinate the wanted processor number. Gets
+        //- Return for every coordinate the wanted processor number.
-        //  passed agglomeration map (from fine to coarse cells) and coarse cell
-        //  location. Can be overridden by decomposers that provide this
-        //  functionality natively.
-        //  See note on weights above.
         virtual labelList decompose
         (
             const polyMesh& mesh,
@@ -308,14 +299,7 @@ public:
             const scalarField& regionWeights
         );
 
-        //- Return for every coordinate the wanted processor number. Explicitly
+        //- Return for every coordinate the wanted processor number.
-        //  provided mesh connectivity.
-        //  The connectivity is equal to mesh.cellCells() except for
-        //  - in parallel the cell numbers are global cell numbers (starting
-        //    from 0 at processor0 and then incrementing all through the
-        //    processors)
-        //  - the connections are across coupled patches
-        //  See note on weights above.
         virtual labelList decompose
         (
             const labelListList& globalCellCells,