STYLE: minor cleanup in decompositionMethod

- reduce cell looping. Avoid initial looping over blockFace.

- make early return (no processor sets, connections, or blocked faces)
  more apparent.

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

@@ -31,6 +31,7 @@ License
 #include "regionSplit.H"
 #include "localPointRegion.H"
 #include "minData.H"
+#include "BitOps.H"
 #include "FaceCellWave.H"
 
 #include "preserveBafflesConstraint.H"
@@ -48,7 +49,8 @@ namespace Foam
 
     // Fallback name when searching for optional coefficients directories
     static const word defaultName("coeffs");
-}
+
+} // End namespace Foam
 
 
 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
@@ -654,15 +656,15 @@ void Foam::decompositionMethod::calcCellCells
     // Check for duplicates connections between cells
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     // Done as postprocessing step since we now have cellCells.
-    label newIndex = 0;
-    labelHashSet nbrCells;
-
 
     if (cellCells.size() == 0)
     {
         return;
     }
 
+    label newIndex = 0;
+    labelHashSet nbrCells;
+
     label startIndex = cellCells.offsets()[0];
 
     forAll(cellCells, celli)
@@ -783,8 +785,8 @@ void Foam::decompositionMethod::calcCellCells
             forAll(pp, i)
             {
                 const label own = agglom[faceOwner[facei]];
-
                 const label globalNei = globalNeighbour[bFacei];
+
                 if
                 (
                    !globalAgglom.isLocal(globalNei)
@@ -839,7 +841,6 @@ void Foam::decompositionMethod::calcCellCells
             forAll(pp, i)
             {
                 const label own = agglom[faceOwner[facei]];
-
                 const label globalNei = globalNeighbour[bFacei];
 
                 if
@@ -848,7 +849,7 @@ void Foam::decompositionMethod::calcCellCells
                  || globalAgglom.toLocal(globalNei) != own
                 )
                 {
-                    label ownIndex = offsets[own] + nFacesPerCell[own]++;
+                    const label ownIndex = offsets[own] + nFacesPerCell[own]++;
                     m[ownIndex] = globalNei;
                     w[ownIndex] = mag(mesh.faceAreas()[facei]);
                 }
@@ -863,15 +864,15 @@ void Foam::decompositionMethod::calcCellCells
     // Check for duplicates connections between cells
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     // Done as postprocessing step since we now have cellCells.
-    label newIndex = 0;
-    labelHashSet nbrCells;
-
 
     if (cellCells.size() == 0)
     {
         return;
     }
 
+    label newIndex = 0;
+    labelHashSet nbrCells;
+
     label startIndex = cellCells.offsets()[0];
 
     forAll(cellCells, celli)
@@ -923,9 +924,9 @@ Foam::labelList Foam::decompositionMethod::decompose
 ) const
 {
     // Any weights specified?
-    label nWeights = returnReduce(cellWeights.size(), sumOp<label>());
+    const bool hasWeights = returnReduce(!cellWeights.empty(), orOp<bool>());
 
-    if (nWeights > 0 && cellWeights.size() != mesh.nCells())
+    if (hasWeights && cellWeights.size() != mesh.nCells())
     {
         FatalErrorInFunction
             << "Number of weights " << cellWeights.size()
@@ -933,115 +934,104 @@ Foam::labelList Foam::decompositionMethod::decompose
             << exit(FatalError);
     }
 
+    // Any faces not blocked?
+    const bool hasUnblocked =
+        returnReduce
+        (
+            (!blockedFace.empty() && !BitOps::all(blockedFace)),
+            orOp<bool>()
+        );
 
-    // Any processor sets?
-    label nProcSets = 0;
-    forAll(specifiedProcessorFaces, setI)
-    {
-        nProcSets += specifiedProcessorFaces[setI].size();
-    }
-    reduce(nProcSets, sumOp<label>());
 
     // Any non-mesh connections?
-    label nConnections = returnReduce
+    const label nConnections = returnReduce
     (
         explicitConnections.size(),
         sumOp<label>()
     );
 
-    // Any faces not blocked?
-    label nUnblocked = 0;
-    forAll(blockedFace, facei)
-    {
-        if (!blockedFace[facei])
-        {
-            nUnblocked++;
-        }
-    }
-    reduce(nUnblocked, sumOp<label>());
 
+    // Any processor sets?
+    label nProcSets = 0;
+    for (const labelList& procset : specifiedProcessorFaces)
+    {
+        nProcSets += procset.size();
+    }
+    reduce(nProcSets, sumOp<label>());
 
 
     // Either do decomposition on cell centres or on agglomeration
 
-    labelList finalDecomp;
+    if (!hasUnblocked && !nConnections && !nProcSets)
-
-
-    if (nProcSets+nConnections+nUnblocked == 0)
     {
         // No constraints, possibly weights
 
-        if (nWeights > 0)
+        return
-        {
-            finalDecomp = decompose
         (
-                mesh,
+            hasWeights
-                mesh.cellCentres(),
+          ? decompose(mesh, mesh.cellCentres(), cellWeights)
-                cellWeights
+          : decompose(mesh, mesh.cellCentres())
         );
     }
-        else
+
-        {
+
-            finalDecomp = decompose(mesh, mesh.cellCentres());
+    // The harder work.
-        }
+    // When we have processor sets, connections, or blocked faces.
-    }
+
-    else
-    {
-        if (debug)
-        {
-            Info<< "Constrained decomposition:" << endl
-                << "    faces with same owner and neighbour processor : "
-                << nUnblocked << endl
-                << "    baffle faces with same owner processor        : "
-                << nConnections << endl
-                << "    faces all on same processor                   : "
-                << nProcSets << endl << endl;
-        }
 
     // Determine local regions, separated by blockedFaces
     regionSplit localRegion(mesh, blockedFace, explicitConnections, false);
 
-
     if (debug)
     {
-            Info<< "Constrained decomposition:" << endl
+        // Only need to count unblocked faces for debugging
+        const label nUnblocked =
+        (
+            hasUnblocked
+          ? returnReduce
+            (
+                label(BitOps::count(blockedFace, false)),
+                sumOp<label>()
+            )
+          : 0
+        );
+
+        Info<< "Constrained decomposition:" << nl
+            << "    faces with same owner and neighbour processor : "
+            << nUnblocked << nl
+            << "    baffle faces with same owner processor        : "
+            << nConnections << nl
+            << "    faces all on same processor                   : "
+            << nProcSets << nl
             << "    split into " << localRegion.nLocalRegions()
             << " regions."
             << endl;
     }
 
-        // Determine region cell centres
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 
-        // This just takes the first cell in the region. Otherwise the problem
+    // Gather region weights and determine region cell centres
-        // is with cyclics - if we'd average the region centre might be
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        // somewhere in the middle of the domain which might not be anywhere
+
-        // near any of the cells.
+    // For the region centre, just take the first cell in the region.
+    // If we average the region centre instead, cyclics could cause
+    // the average domain centre to be outside of domain.
+
+    scalarField regionWeights(localRegion.nLocalRegions(), 0.0);
 
     pointField regionCentres(localRegion.nLocalRegions(), point::max);
 
-        forAll(localRegion, celli)
+    if (hasWeights)
-        {
-            const label regioni = localRegion[celli];
-
-            if (regionCentres[regioni] == point::max)
-            {
-                regionCentres[regioni] = mesh.cellCentres()[celli];
-            }
-        }
-
-        // Do decomposition on agglomeration
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
-        scalarField regionWeights(localRegion.nLocalRegions(), 0);
-
-        if (nWeights > 0)
     {
         forAll(localRegion, celli)
         {
             const label regioni = localRegion[celli];
 
             regionWeights[regioni] += cellWeights[celli];
+
+            if (regionCentres[regioni] == point::max)
+            {
+                regionCentres[regioni] = mesh.cellCentres()[celli];
+            }
         }
     }
     else
@@ -1051,10 +1041,19 @@ Foam::labelList Foam::decompositionMethod::decompose
             const label regioni = localRegion[celli];
 
             regionWeights[regioni] += 1.0;
+
+            if (regionCentres[regioni] == point::max)
+            {
+                regionCentres[regioni] = mesh.cellCentres()[celli];
+            }
         }
     }
 
-        finalDecomp = decompose
+    // Do decomposition on agglomeration
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    labelList finalDecomp =
+        decompose
         (
             mesh,
             localRegion,
@@ -1063,20 +1062,17 @@ Foam::labelList Foam::decompositionMethod::decompose
         );
 
 
-
+    // Apply explicitConnections since decompose did not know about them
-        // Implement the explicitConnections since above decompose
+    for (const labelPair& baffle : explicitConnections)
-        // does not know about them
-        forAll(explicitConnections, connectioni)
     {
-            const labelPair& baffle = explicitConnections[connectioni];
         const label f0 = baffle.first();
         const label f1 = baffle.second();
 
         if (!blockedFace[f0] && !blockedFace[f1])
         {
             // Note: what if internal faces and owner and neighbour on
-                // different processor? So for now just push owner side
+            // different processor?
-                // proc
+            // So for now just push owner side proc
 
             const label proci = finalDecomp[mesh.faceOwner()[f0]];
 
@@ -1102,6 +1098,9 @@ Foam::labelList Foam::decompositionMethod::decompose
     // separately since not handled by local regionSplit. We need to
     // walk now across coupled faces and make sure to move a whole
     // global region across
+
+    // This additionally consolidates/compacts the regions numbers globally,
+    // since that was skipped in the previous regionSplit.
     if (Pstream::parRun())
     {
         // Re-do regionSplit
@@ -1121,7 +1120,7 @@ Foam::labelList Foam::decompositionMethod::decompose
             }
             else
             {
-                    nUnblocked++;
+                ++nUnblocked;
             }
         }
 
@@ -1177,15 +1176,14 @@ Foam::labelList Foam::decompositionMethod::decompose
     //
     // Note that reworking the cellToProc might make the decomposition
     // unbalanced.
-        forAll(specifiedProcessorFaces, setI)
+    forAll(specifiedProcessorFaces, seti)
     {
-            const labelList& set = specifiedProcessorFaces[setI];
+        const labelList& set = specifiedProcessorFaces[seti];
 
-            label proci = specifiedProcessor[setI];
+        label proci = specifiedProcessor[seti];
         if (proci == -1)
         {
-                // If no processor specified use the one from the
+            // If no processor specified - use the one from the 0th element
-                // 0th element
             if (set.size())
             {
                 proci = finalDecomp[mesh.faceOwner()[set[0]]];
@@ -1197,18 +1195,18 @@ Foam::labelList Foam::decompositionMethod::decompose
             }
         }
 
-            forAll(set, fI)
+        for (const label facei : set)
         {
-                const face& f = mesh.faces()[set[fI]];
+            const face& f = mesh.faces()[facei];
-                forAll(f, fp)
+            for (const label pointi : f)
             {
-                    const labelList& pFaces = mesh.pointFaces()[f[fp]];
+                const labelList& pFaces = mesh.pointFaces()[pointi];
-                    for (const label facei : pFaces)
+                for (const label pFacei : pFaces)
                 {
-                        finalDecomp[mesh.faceOwner()[facei]] = proci;
+                    finalDecomp[mesh.faceOwner()[pFacei]] = proci;
-                        if (mesh.isInternalFace(facei))
+                    if (mesh.isInternalFace(pFacei))
                     {
-                            finalDecomp[mesh.faceNeighbour()[facei]] = proci;
+                        finalDecomp[mesh.faceNeighbour()[pFacei]] = proci;
                     }
                 }
             }
@@ -1228,14 +1226,15 @@ Foam::labelList Foam::decompositionMethod::decompose
             {
                 forAll(pp, i)
                 {
-                        label facei = pp.start()+i;
+                    const label facei = pp.start()+i;
-                        label own = mesh.faceOwner()[facei];
+                    const label own = mesh.faceOwner()[facei];
-                        label bFacei = facei-mesh.nInternalFaces();
+                    const label bFacei = facei-mesh.nInternalFaces();
 
                     if (!blockedFace[facei])
                     {
-                            label ownProc = finalDecomp[own];
+                        const label ownProc = finalDecomp[own];
-                            label nbrProc = nbrDecomp[bFacei];
+                        const label nbrProc = nbrDecomp[bFacei];
+
                         if (ownProc != nbrProc)
                         {
                             FatalErrorInFunction
@@ -1251,7 +1250,6 @@ Foam::labelList Foam::decompositionMethod::decompose
             }
         }
     }
-    }
 
     return finalDecomp;
 }