BUG: #1479: faceZone inside cellZone

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

@@ -2353,7 +2353,21 @@ Foam::labelList Foam::meshRefinement::freeStandingBaffleFaces
     const labelList& faceNeighbour = mesh_.faceNeighbour();
 
 
-    DynamicList<label> faceLabels(mesh_.nFaces()/20);
+    // We want to pick up the faces to orient. These faces come in
+    // two variants:
+    // - faces originating from stand-alone faceZones
+    //   (these will most likely have no cellZone on either side so
+    //    ownZone and neiZone both -1)
+    // - sticky-up faces originating from a 'bulge' in a outside of
+    //   a cellZone. These will have the same cellZone on either side.
+    //   How to orient these is not really clearly defined so do them
+    //   same as stand-alone faceZone faces for now. (Normally these will
+    //   already have been removed by the 'allowFreeStandingZoneFaces=false'
+    //   default setting)
+
+    // Note that argument neiCellZone will have -1 on uncoupled boundaries.
+
+    DynamicList<label> faceLabels(mesh_.nFaces()/100);
 
     for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
     {
@@ -2362,7 +2376,7 @@ Foam::labelList Foam::meshRefinement::freeStandingBaffleFaces
             // Free standing baffle?
             label ownZone = cellToZone[faceOwner[faceI]];
             label neiZone = cellToZone[faceNeighbour[faceI]];
-            if (max(ownZone, neiZone) == -1)
+            if (ownZone == neiZone)
             {
                 faceLabels.append(faceI);
             }
@@ -2380,7 +2394,7 @@ Foam::labelList Foam::meshRefinement::freeStandingBaffleFaces
                 // Free standing baffle?
                 label ownZone = cellToZone[faceOwner[faceI]];
                 label neiZone = neiCellZone[faceI-mesh_.nInternalFaces()];
-                if (max(ownZone, neiZone) == -1)
+                if (ownZone == neiZone)
                 {
                     faceLabels.append(faceI);
                 }
@@ -2853,12 +2867,15 @@ void Foam::meshRefinement::zonify
         if (faceZoneI != -1)
         {
             // Orient face zone to have slave cells in min cell zone.
+            // Note: logic to use flipMap should be consistent with logic
+            //       to pick up the freeStandingBaffleFaces!
+
             label ownZone = cellToZone[faceOwner[faceI]];
             label neiZone = cellToZone[faceNeighbour[faceI]];
 
             bool flip;
 
-            if (ownZone == -1 && neiZone == -1)
+            if (ownZone == neiZone)
             {
                 // free-standing face. Use geometrically derived orientation
                 flip = meshFlipMap[faceI];
@@ -2911,17 +2928,11 @@ void Foam::meshRefinement::zonify
 
                 bool flip;
 
-                if (ownZone == -1 && neiZone == -1)
+                if (ownZone == neiZone)
                 {
                     // free-standing face. Use geometrically derived orientation
                     flip = meshFlipMap[faceI];
                 }
-                else if (ownZone == neiZone)
-                {
-                    // Free-standing zone face or coupled boundary. Keep
-                    // master face unflipped.
-                    flip = !isMasterFace[faceI];
-                }
                 else
                 {
                     flip =
@@ -4383,22 +4394,23 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify
 
 
     // Get coupled neighbour cellZone. Set to -1 on non-coupled patches.
-    labelList neiCellZone(mesh_.nFaces()-mesh_.nInternalFaces(), -1);
+    labelList neiCellZone;
+    syncTools::swapBoundaryCellList(mesh_, cellToZone, neiCellZone);
     forAll(patches, patchI)
     {
         const polyPatch& pp = patches[patchI];
 
-        if (pp.coupled())
+        if (!pp.coupled())
         {
+            label bFaceI = pp.start()-mesh_.nInternalFaces();
             forAll(pp, i)
             {
-                label faceI = pp.start()+i;
-                neiCellZone[faceI-mesh_.nInternalFaces()] =
-                    cellToZone[mesh_.faceOwner()[faceI]];
+                neiCellZone[bFaceI++] = -1;
             }
         }
     }
-    syncTools::swapBoundaryFaceList(mesh_, neiCellZone);
+
+
 
     // Get per face whether is it master (of a coupled set of faces)
     const PackedBoolList isMasterFace(syncTools::getMasterFaces(mesh_));
@@ -4438,6 +4450,13 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify
             Info<< "Detected " << nFreeStanding << " free-standing zone faces"
                 << endl;
 
+            if (debug)
+            {
+                OBJstream str(mesh_.time().path()/"freeStanding.obj");
+                str.write(patch.localFaces(), patch.localPoints(), false);
+            }
+
+
             // Detect non-manifold edges
             labelList nMasterFacesPerEdge;
             calcPatchNumMasterFaces(isMasterFace, patch, nMasterFacesPerEdge);