ENH: Do a weighted balance before refinement and layer addition.

Have weight=7 on cells to be refined and balance. Then after refinement
there will be perfect balancing. Similar for layer addition.

applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C

@@ -434,6 +434,19 @@ int main(int argc, char *argv[])
         // Layer addition parameters
         layerParameters layerParams(layerDict, mesh.boundaryMesh());
 
+        //!!! Temporary hack to get access to maxLocalCells
+        bool preBalance;
+        {
+            refinementParameters refineParams(refineDict);
+
+            preBalance = returnReduce
+            (
+                (mesh.nCells() >= refineParams.maxLocalCells()),
+                orOp<bool>()
+            );
+        }
+
+
         if (!overwrite)
         {
             const_cast<Time&>(mesh.time())++;
@@ -444,6 +457,7 @@ int main(int argc, char *argv[])
             layerDict,
             motionDict,
             layerParams,
+            preBalance,
             decomposer,
             distributor
         );

applications/utilities/mesh/generation/snappyHexMesh/snappyHexMeshDict

@@ -10,7 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       dictionary;
-    object      autoHexMeshDict;
+    object      snappyHexMeshDict;
 }
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -40,7 +40,7 @@ geometry
     {
         type triSurfaceMesh;
 
-        //tolerance   1E-6;   // optional:non-default tolerance on intersections
+        //tolerance   1E-5;   // optional:non-default tolerance on intersections
         //maxTreeDepth 10;    // optional:depth of octree. Decrease only in case
                               // of memory limitations.
 
@@ -71,9 +71,9 @@ castellatedMeshControls
     // Refinement parameters
     // ~~~~~~~~~~~~~~~~~~~~~
 
-    // While refining maximum number of cells per processor. This is basically
+    // If local number of cells is >= maxLocalCells on any processor
-    // the number of cells that fit on a processor. If you choose this too small
+    // switches from from refinement followed by balancing
-    // it will do just more refinement iterations to obtain a similar mesh.
+    // (current method) to (weighted) balancing before refinement.
     maxLocalCells 1000000;
 
     // Overall cell limit (approximately). Refinement will stop immediately
@@ -89,6 +89,13 @@ castellatedMeshControls
     // (unless the number of cells to refine is 0)
     minRefinementCells 0;
 
+    // Allow a certain level of imbalance during refining
+    // (since balancing is quite expensive)
+    // Expressed as fraction of perfect balance (= overall number of cells /
+    // nProcs). 0=balance always.
+    maxLoadUnbalance 0.10;
+
+
     // Number of buffer layers between different levels.
     // 1 means normal 2:1 refinement restriction, larger means slower
     // refinement.

src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.C

@@ -539,6 +539,7 @@ void Foam::autoHexMeshDriver::doMesh()
             shrinkDict,
             motionDict,
             layerParams,
+            true,                   // pre-balance
             decomposerPtr_(),
             distributorPtr_()
         );

src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.C

@@ -2497,18 +2497,13 @@ void Foam::autoLayerDriver::mergePatchFacesUndo
         << "---------------------------" << nl
         << "    - which are on the same patch" << nl
         << "    - which make an angle < " << layerParams.featureAngle()
-        << "- which are on the same patch" << nl
-        << "- which make an angle < " << layerParams.featureAngle()
         << " degrees"
         << nl
         << "      (cos:" << minCos << ')' << nl
         << "    - as long as the resulting face doesn't become concave"
-        << "  (cos:" << minCos << ')' << nl
-        << "- as long as the resulting face doesn't become concave"
         << " by more than "
         << layerParams.concaveAngle() << " degrees" << nl
         << "      (0=straight, 180=fully concave)" << nl
-        << "  (0=straight, 180=fully concave)" << nl
         << endl;
 
     label nChanged = mergePatchFacesUndo(minCos, concaveCos, motionDict);
@@ -3241,6 +3236,7 @@ void Foam::autoLayerDriver::doLayers
     const dictionary& shrinkDict,
     const dictionary& motionDict,
     const layerParameters& layerParams,
+    const bool preBalance,
     decompositionMethod& decomposer,
     fvMeshDistribute& distributor
 )
@@ -3299,6 +3295,7 @@ void Foam::autoLayerDriver::doLayers
 
         // Balance
         if (Pstream::parRun())
+        if (Pstream::parRun() && preBalance)
         {
             Info<< nl
                 << "Doing initial balancing" << nl

src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.H

@@ -535,6 +535,7 @@ public:
                 const dictionary& shrinkDict,
                 const dictionary& motionDict,
                 const layerParameters& layerParams,
+                const bool preBalance,              // balance before adding?
                 decompositionMethod& decomposer,
                 fvMeshDistribute& distributor
             );

src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoRefineDriver.C

@@ -205,14 +205,36 @@ Foam::label Foam::autoRefineDriver::featureEdgeRefine
                 const_cast<Time&>(mesh.time())++;
             }
 
-            meshRefiner_.refineAndBalance
+
+            if
             (
-                "feature refinement iteration " + name(iter),
+                returnReduce
-                decomposer_,
+                (
-                distributor_,
+                    (mesh.nCells() >= refineParams.maxLocalCells()),
-                cellsToRefine,
+                    orOp<bool>()
-                refineParams.maxLoadUnbalance()
+                )
-            );
+            )
+            {
+                meshRefiner_.balanceAndRefine
+                (
+                    "feature refinement iteration " + name(iter),
+                    decomposer_,
+                    distributor_,
+                    cellsToRefine,
+                    refineParams.maxLoadUnbalance()
+                );
+            }
+            else
+            {
+                meshRefiner_.refineAndBalance
+                (
+                    "feature refinement iteration " + name(iter),
+                    decomposer_,
+                    distributor_,
+                    cellsToRefine,
+                    refineParams.maxLoadUnbalance()
+                );
+            }
         }
     }
     return iter;
@@ -306,14 +328,36 @@ Foam::label Foam::autoRefineDriver::surfaceOnlyRefine
             const_cast<Time&>(mesh.time())++;
         }
 
-        meshRefiner_.refineAndBalance
+
+        if
         (
-            "surface refinement iteration " + name(iter),
+            returnReduce
-            decomposer_,
+            (
-            distributor_,
+                (mesh.nCells() >= refineParams.maxLocalCells()),
-            cellsToRefine,
+                orOp<bool>()
-            refineParams.maxLoadUnbalance()
+            )
-        );
+        )
+        {
+            meshRefiner_.balanceAndRefine
+            (
+                "surface refinement iteration " + name(iter),
+                decomposer_,
+                distributor_,
+                cellsToRefine,
+                refineParams.maxLoadUnbalance()
+            );
+        }
+        else
+        {
+            meshRefiner_.refineAndBalance
+            (
+                "surface refinement iteration " + name(iter),
+                decomposer_,
+                distributor_,
+                cellsToRefine,
+                refineParams.maxLoadUnbalance()
+            );
+        }
     }
     return iter;
 }
@@ -491,14 +535,35 @@ Foam::label Foam::autoRefineDriver::shellRefine
             const_cast<Time&>(mesh.time())++;
         }
 
-        meshRefiner_.refineAndBalance
+        if
         (
-            "shell refinement iteration " + name(iter),
+            returnReduce
-            decomposer_,
+            (
-            distributor_,
+                (mesh.nCells() >= refineParams.maxLocalCells()),
-            cellsToRefine,
+                orOp<bool>()
-            refineParams.maxLoadUnbalance()
+            )
-        );
+        )
+        {
+            meshRefiner_.balanceAndRefine
+            (
+                "shell refinement iteration " + name(iter),
+                decomposer_,
+                distributor_,
+                cellsToRefine,
+                refineParams.maxLoadUnbalance()
+            );
+        }
+        else
+        {
+            meshRefiner_.refineAndBalance
+            (
+                "shell refinement iteration " + name(iter),
+                decomposer_,
+                distributor_,
+                cellsToRefine,
+                refineParams.maxLoadUnbalance()
+            );
+        }
     }
     meshRefiner_.userFaceData().clear();
 

src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinement.H

@@ -656,6 +656,16 @@ public:
                 const scalar maxLoadUnbalance
             );
 
+            //- Balance before refining some cells
+            autoPtr<mapDistributePolyMesh> balanceAndRefine
+            (
+                const string& msg,
+                decompositionMethod& decomposer,
+                fvMeshDistribute& distributor,
+                const labelList& cellsToRefine,
+                const scalar maxLoadUnbalance
+            );
+
 
         // Baffle handling
 

src/mesh/autoMesh/autoHexMesh/meshRefinement/meshRefinementRefine.C

@@ -219,7 +219,6 @@ Foam::labelList Foam::meshRefinement::getChangedFaces
             << endl;
 
         faceSet changedFacesSet(mesh, "changedFaces", changedFaces);
-        changedFacesSet.instance() = mesh.time().timeName();
         Pout<< "getChangedFaces : Writing " << changedFaces.size()
             << " changed faces to faceSet " << changedFacesSet.name()
             << endl;
@@ -1246,90 +1245,110 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::refine
 }
 
 
-//// Do refinement of consistent set of cells followed by truncation and
+// Do refinement of consistent set of cells followed by truncation and
-//// load balancing.
+// load balancing.
-//Foam::autoPtr<Foam::mapDistributePolyMesh>
+Foam::autoPtr<Foam::mapDistributePolyMesh>
-//Foam::meshRefinement::refineAndBalance
+Foam::meshRefinement::refineAndBalance
-//(
+(
-//    const string& msg,
+    const string& msg,
-//    decompositionMethod& decomposer,
+    decompositionMethod& decomposer,
-//    fvMeshDistribute& distributor,
+    fvMeshDistribute& distributor,
-//    const labelList& cellsToRefine
+    const labelList& cellsToRefine,
-//)
+    const scalar maxLoadUnbalance
-//{
+)
-//    // Do all refinement
+{
-//    refine(cellsToRefine);
+    // Do all refinement
-//
+    refine(cellsToRefine);
-//    if (debug)
+
-//    {
+    if (debug)
-//        Pout<< "Writing refined but unbalanced " << msg
+    {
-//            << " mesh to time " << timeName() << endl;
+        Pout<< "Writing refined but unbalanced " << msg
-//        write
+            << " mesh to time " << timeName() << endl;
-//        (
+        write
-//            debug,
+        (
-//            mesh_.time().path()
+            debug,
-//           /timeName()
+            mesh_.time().path()
-//        );
+           /timeName()
-//        Pout<< "Dumped debug data in = "
+        );
-//            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+        Pout<< "Dumped debug data in = "
-//
+            << mesh_.time().cpuTimeIncrement() << " s" << endl;
-//        // test all is still synced across proc patches
+
-//        checkData();
+        // test all is still synced across proc patches
-//    }
+        checkData();
-//
+    }
-//    Info<< "Refined mesh in = "
+
-//        << mesh_.time().cpuTimeIncrement() << " s" << endl;
+    Info<< "Refined mesh in = "
-//    printMeshInfo(debug, "After refinement " + msg);
+        << mesh_.time().cpuTimeIncrement() << " s" << endl;
-//
+    printMeshInfo(debug, "After refinement " + msg);
-//
+
-//    // Load balancing
+
-//    // ~~~~~~~~~~~~~~
+    // Load balancing
-//
+    // ~~~~~~~~~~~~~~
-//    autoPtr<mapDistributePolyMesh> distMap;
+
-//
+    autoPtr<mapDistributePolyMesh> distMap;
-//    if (Pstream::nProcs() > 1)
+
-//    {
+    if (Pstream::nProcs() > 1)
-//        scalarField cellWeights(mesh_.nCells(), 1);
+    {
-//
+        scalar nIdealCells =
-//        distMap = balance
+            mesh_.globalData().nTotalCells()
-//        (
+          / Pstream::nProcs();
-//            false,  //keepZoneFaces
+
-//            false,  //keepBaffles
+        scalar unbalance = returnReduce
-//            cellWeights,
+        (
-//            decomposer,
+            mag(1.0-mesh_.nCells()/nIdealCells),
-//            distributor
+            maxOp<scalar>()
-//        );
+        );
-//
+
-//        Info<< "Balanced mesh in = "
+        if (unbalance <= maxLoadUnbalance)
-//            << mesh_.time().cpuTimeIncrement() << " s" << endl;
+        {
-//
+            Info<< "Skipping balancing since max unbalance " << unbalance
-//        printMeshInfo(debug, "After balancing " + msg);
+                << " is less than allowable " << maxLoadUnbalance
-//
+                << endl;
-//
+        }
-//        if (debug)
+        else
-//        {
+        {
-//            Pout<< "Writing balanced " << msg
+            scalarField cellWeights(mesh_.nCells(), 1);
-//                << " mesh to time " << timeName() << endl;
+
-//            write
+            distMap = balance
-//            (
+            (
-//                debug,
+                false,  //keepZoneFaces
-//                mesh_.time().path()/timeName()
+                false,  //keepBaffles
-//            );
+                cellWeights,
-//            Pout<< "Dumped debug data in = "
+                decomposer,
-//                << mesh_.time().cpuTimeIncrement() << " s" << endl;
+                distributor
-//
+            );
-//            // test all is still synced across proc patches
+
-//            checkData();
+            Info<< "Balanced mesh in = "
-//        }
+                << mesh_.time().cpuTimeIncrement() << " s" << endl;
-//    }
+
-//
+            printMeshInfo(debug, "After balancing " + msg);
-//    return distMap;
+
-//}
+
+            if (debug)
+            {
+                Pout<< "Writing balanced " << msg
+                    << " mesh to time " << timeName() << endl;
+                write
+                (
+                    debug,
+                    mesh_.time().path()/timeName()
+                );
+                Pout<< "Dumped debug data in = "
+                    << mesh_.time().cpuTimeIncrement() << " s" << endl;
+
+                // test all is still synced across proc patches
+                checkData();
+            }
+        }
+    }
+
+    return distMap;
+}
 
 
 // Do load balancing followed by refinement of consistent set of cells.
 Foam::autoPtr<Foam::mapDistributePolyMesh>
-Foam::meshRefinement::refineAndBalance
+Foam::meshRefinement::balanceAndRefine
 (
     const string& msg,
     decompositionMethod& decomposer,
@@ -1386,8 +1405,8 @@ Foam::meshRefinement::refineAndBalance
         if (unbalance <= maxLoadUnbalance)
         {
             Info<< "Skipping balancing since max unbalance " << unbalance
-                << " in = "
+                << " is less than allowable " << maxLoadUnbalance
-                << mesh_.time().cpuTimeIncrement() << " s" << endl;
+                << endl;
         }
         else
         {

tutorials/heatTransfer/buoyantBoussinesqSimpleFoam/iglooWithFridges/system/snappyHexMeshDict

@@ -116,9 +116,9 @@ castellatedMeshControls
     // Refinement parameters
     // ~~~~~~~~~~~~~~~~~~~~~
 
-    // While refining maximum number of cells per processor. This is basically
+    // If local number of cells is >= maxLocalCells on any processor
-    // the number of cells that fit on a processor. If you choose this too small
+    // switches from from refinement followed by balancing
-    // it will do just more refinement iterations to obtain a similar mesh.
+    // (current method) to (weighted) balancing before refinement.
     maxLocalCells 1000000;
 
     // Overall cell limit (approximately). Refinement will stop immediately
@@ -255,6 +255,8 @@ snapControls
 // Settings for the layer addition.
 addLayersControls
 {
+    // Are the thickness parameters below relative to the undistorted
+    // size of the refined cell outside layer (true) or absolute sizes (false).
     relativeSizes true;
 
     // Per final patch (so not geometry!) the layer information
@@ -277,11 +279,14 @@ addLayersControls
     //  is the
     //  thickness of the layer furthest away from the wall.
     //  Relative to undistorted size of cell outside layer.
+    //  is the thickness of the layer furthest away from the wall.
+    //  See relativeSizes parameter.
     finalLayerThickness 0.5;
 
     //- Minimum thickness of cell layer. If for any reason layer
     //  cannot be above minThickness do not add layer.
     //  Relative to undistorted size of cell outside layer.
+    //  See relativeSizes parameter.
     minThickness 0.25;
 
     //- If points get not extruded do nGrow layers of connected faces that are
@@ -323,7 +328,10 @@ addLayersControls
     // Create buffer region for new layer terminations
     nBufferCellsNoExtrude 0;
 
-    // Overall max number of layer addition iterations
+
+    // Overall max number of layer addition iterations. The mesher will exit
+    // if it reaches this number of iterations; possibly with an illegal
+    // mesh.
     nLayerIter 50;
 }
 
@@ -349,6 +357,7 @@ meshQualityControls
     minFlatness 0.5;
 
     //- Minimum pyramid volume. Is absolute volume of cell pyramid.
+    //  Set to a sensible fraction of the smallest cell volume expected.
     //  Set to very negative number (e.g. -1E30) to disable.
     minVol 1e-13;
 

tutorials/heatTransfer/chtMultiRegionFoam/snappyMultiRegionHeater/system/snappyHexMeshDict

@@ -210,6 +210,12 @@ castellatedMeshControls
     // NOTE: This point should never be on a face, always inside a cell, even
     // after refinement.
     locationInMesh (0.01 0.01 0.01);
+
+
+    // Whether any faceZones (as specified in the refinementSurfaces)
+    // are only on the boundary of corresponding cellZones or also allow
+    // free-standing zone faces. Not used if there are no faceZones.
+    allowFreeStandingZoneFaces false;
 }
 
 

tutorials/incompressible/simpleFoam/motorBike/system/snappyHexMeshDict

@@ -51,9 +51,9 @@ castellatedMeshControls
     // Refinement parameters
     // ~~~~~~~~~~~~~~~~~~~~~
 
-    // While refining maximum number of cells per processor. This is basically
+    // If local number of cells is >= maxLocalCells on any processor
-    // the number of cells that fit on a processor. If you choose this too small
+    // switches from from refinement followed by balancing
-    // it will do just more refinement iterations to obtain a similar mesh.
+    // (current method) to (weighted) balancing before refinement.
     maxLocalCells 1000000;
 
     // Overall cell limit (approximately). Refinement will stop immediately
@@ -69,6 +69,13 @@ castellatedMeshControls
     // (unless the number of cells to refine is 0)
     minRefinementCells 10;
 
+    // Allow a certain level of imbalance during refining
+    // (since balancing is quite expensive)
+    // Expressed as fraction of perfect balance (= overall number of cells /
+    // nProcs). 0=balance always.
+    maxLoadUnbalance 0.10;
+
+
     // Number of buffer layers between different levels.
     // 1 means normal 2:1 refinement restriction, larger means slower
     // refinement.
@@ -180,6 +187,8 @@ snapControls
 // Settings for the layer addition.
 addLayersControls
 {
+    // Are the thickness parameters below relative to the undistorted
+    // size of the refined cell outside layer (true) or absolute sizes (false).
     relativeSizes true;
 
     // Per final patch (so not geometry!) the layer information
@@ -466,6 +475,8 @@ addLayersControls
     //  is the
     //  thickness of the layer furthest away from the wall.
     //  Relative to undistorted size of cell outside layer.
+    //  is the thickness of the layer furthest away from the wall.
+    //  See relativeSizes parameter.
     finalLayerThickness 0.3;
 
     //- Minimum thickness of cell layer. If for any reason layer
@@ -539,6 +550,7 @@ meshQualityControls
     minFlatness 0.5;
 
     //- Minimum pyramid volume. Is absolute volume of cell pyramid.
+    //  Set to a sensible fraction of the smallest cell volume expected.
     //  Set to very negative number (e.g. -1E30) to disable.
     minVol 1e-13;