Merge branch 'master' of ssh://noisy/home/noisy2/OpenFOAM/OpenFOAM-dev

applications/utilities/mesh/conversion/netgenNeutralToFoam/netgenNeutralToFoam.C

@@ -89,13 +89,11 @@ using namespace Foam;
 int main(int argc, char *argv[])
 {
     argList::validArgs.append("Neutral file");
-    argList::validOptions.insert("overwrite", "");
 
 #   include "setRootCase.H"
 #   include "createTime.H"
 
     fileName neuFile(args.additionalArgs()[0]);
-    bool overwrite = args.options().found("overwrite");
 
 
     IFstream str(neuFile);
@@ -300,11 +298,6 @@ int main(int argc, char *argv[])
     }
 
 
-    if (!overwrite)
-    {
-        runTime++;
-    }
-
     polyMesh mesh
     (
         IOobject

applications/utilities/mesh/conversion/tetgenToFoam/tetgenToFoam.C

@@ -54,8 +54,9 @@ Description
 NOTE:
 - for some reason boundary faces point inwards. I just reverse them
 always. Might use some geometric check instead.
-- marked faces might not actually be boundary faces of mesh. This is not handled
-and you'll have to run without face file (-noFaceFile option)
+- marked faces might not actually be boundary faces of mesh.
+This is hopefully handled now by first constructing without boundaries
+and then reconstructing with boundary faces.
 
 \*---------------------------------------------------------------------------*/
 
@@ -69,20 +70,40 @@ using namespace Foam;
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+// Find label of face.
+label findFace(const primitiveMesh& mesh, const face& f)
+{
+    const labelList& pFaces = mesh.pointFaces()[f[0]];
+
+    forAll(pFaces, i)
+    {
+        label faceI = pFaces[i];
+
+        if (mesh.faces()[faceI] == f)
+        {
+            return faceI;
+        }
+    }
+
+    FatalErrorIn("findFace(const primitiveMesh&, const face&)")
+        << "Cannot find face " << f << " in mesh." << abort(FatalError);
+
+    return -1;
+}
+
+
 // Main program:
 
 int main(int argc, char *argv[])
 {
     argList::validArgs.append("file prefix");
     argList::validOptions.insert("noFaceFile", "");
-    argList::validOptions.insert("overwrite", "");
 
 #   include "setRootCase.H"
 #   include "createTime.H"
 
 
     bool readFaceFile = !args.options().found("noFaceFile");
-    bool overwrite = args.options().found("overwrite");
 
     fileName prefix(args.additionalArgs()[0]);
 
@@ -287,16 +308,44 @@ int main(int argc, char *argv[])
     }
 
 
-    label nPatches = 0;
+    //
+    // Construct mesh with default boundary only
+    //
+
+    autoPtr<polyMesh> meshPtr
+    (
+        new polyMesh
+        (
+            IOobject
+            (
+                polyMesh::defaultRegion,
+                runTime.constant(),
+                runTime
+            ),
+            points,
+            cells,
+            faceListList(0),
+            wordList(0),    //boundaryPatchNames
+            wordList(0),    //boundaryPatchTypes
+            "defaultFaces",
+            polyPatch::typeName,
+            wordList(0)
+        )
+    );
+    const polyMesh& mesh = meshPtr;
 
-    // List of Foam vertices per boundary face
-    faceList boundaryFaces;
 
-    // For each boundary faces the Foam patchID
-    labelList boundaryPatch;
 
     if (readFaceFile)
     {
+        label nPatches = 0;
+
+        // List of Foam vertices per boundary face
+        faceList boundaryFaces;
+
+        // For each boundary faces the Foam patchID
+        labelList boundaryPatch;
+
         //
         // read boundary faces
         //
@@ -366,48 +415,59 @@ int main(int argc, char *argv[])
                     f[2-i] = nodeToPoint[nodeI];
                 }
 
-                boundaryFaces[faceI] = f;
 
-                if (nFaceAttr > 0)
+                if (findFace(mesh, f) >= mesh.nInternalFaces())
                 {
-                    // First attribute is the region number
-                    faceLine >> region;
+                    boundaryFaces[faceI] = f;
 
+                    if (nFaceAttr > 0)
+                    {
+                        // First attribute is the region number
+                        faceLine >> region;
 
-                    // Get Foam patchID and update region->patch table.
-                    label patchI = 0;
 
-                    Map<label>::iterator patchFind = regionToPatch.find(region);
+                        // Get Foam patchID and update region->patch table.
+                        label patchI = 0;
 
-                    if (patchFind == regionToPatch.end())
-                    {
-                        patchI = nPatches;
+                        Map<label>::iterator patchFind =
+                            regionToPatch.find(region);
 
-                        Info<< "Mapping tetgen region " << region
-                            << " to Foam patch "
-                            << patchI << endl;
+                        if (patchFind == regionToPatch.end())
+                        {
+                            patchI = nPatches;
 
-                        regionToPatch.insert(region, nPatches++);
-                    }
-                    else
-                    {
-                        patchI = patchFind();
-                    }
+                            Info<< "Mapping tetgen region " << region
+                                << " to Foam patch "
+                                << patchI << endl;
 
-                    boundaryPatch[faceI] = patchI;
+                            regionToPatch.insert(region, nPatches++);
+                        }
+                        else
+                        {
+                            patchI = patchFind();
+                        }
 
-                    // Skip remaining attributes
-                    for (label i = 1; i < nFaceAttr; i++)
-                    {
-                        faceLine >> dummy;
+                        boundaryPatch[faceI] = patchI;
+
+                        // Skip remaining attributes
+                        for (label i = 1; i < nFaceAttr; i++)
+                        {
+                            faceLine >> dummy;
+                        }
                     }
-                }
 
-                faceI++;
+                    faceI++;
+                }
             }
         }
 
-        // Print region to patch mapping
+   
+        // Trim
+        boundaryFaces.setSize(faceI);
+        boundaryPatch.setSize(faceI);
+
+
+         // Print region to patch mapping
         Info<< "Regions:" << endl;
 
         for
@@ -421,28 +481,23 @@ int main(int argc, char *argv[])
                 << iter() << endl;
         }
         Info<< endl;
-    }
 
 
-    // Storage for boundary faces
+        // Storage for boundary faces
+        faceListList patchFaces(nPatches);
+        wordList patchNames(nPatches);
 
-    faceListList patchFaces(nPatches);
-
-    wordList patchNames(nPatches);
-
-    forAll(patchNames, patchI)
-    {
-        patchNames[patchI] = word("patch") + name(patchI);
-    }
+        forAll(patchNames, patchI)
+        {
+            patchNames[patchI] = word("patch") + name(patchI);
+        }
 
-    wordList patchTypes(nPatches, polyPatch::typeName);
-    word defaultFacesName = "defaultFaces";
-    word defaultFacesType = polyPatch::typeName;
-    wordList patchPhysicalTypes(nPatches, polyPatch::typeName);
+        wordList patchTypes(nPatches, polyPatch::typeName);
+        word defaultFacesName = "defaultFaces";
+        word defaultFacesType = polyPatch::typeName;
+        wordList patchPhysicalTypes(nPatches, polyPatch::typeName);
 
 
-    if (readFaceFile)
-    {
         // Sort boundaryFaces by patch using boundaryPatch.
         List<DynamicList<face> > allPatchFaces(nPatches);
 
@@ -464,34 +519,34 @@ int main(int argc, char *argv[])
         }
 
         Info<< endl;
-    }
 
-    if (!overwrite)
-    {
-        runTime++;
-    }
 
-    polyMesh mesh
-    (
-        IOobject
+        meshPtr.reset
         (
-            polyMesh::defaultRegion,
-            runTime.constant(),
-            runTime
-        ),
-        points,
-        cells,
-        patchFaces,
-        patchNames,
-        patchTypes,
-        defaultFacesName,
-        defaultFacesType,
-        patchPhysicalTypes
-    );
+            new polyMesh
+            (
+                IOobject
+                (
+                    polyMesh::defaultRegion,
+                    runTime.constant(),
+                    runTime
+                ),
+                points,
+                cells,
+                patchFaces,
+                patchNames,
+                patchTypes,
+                defaultFacesName,
+                defaultFacesType,
+                patchPhysicalTypes
+            )
+        );
+    }
+
 
     Info<< "Writing mesh to " << runTime.constant() << endl << endl;
 
-    mesh.write();
+    meshPtr().write();
 
 
     Info<< "End\n" << endl;

src/OpenFOAM/containers/HashTables/HashTable/HashTable.C

@@ -485,6 +485,43 @@ void HashTable<T, Key, Hash>::operator=(const HashTable<T, Key, Hash>& ht)
 }
 
 
+template<class T, class Key, class Hash>
+bool HashTable<T, Key, Hash>::operator==(const HashTable<T, Key, Hash>& ht)
+ const
+{
+    // Are all my elements in ht?
+    for (const_iterator iter = begin(); iter != end(); ++iter)
+    {
+        const_iterator fnd = ht.find(iter.key());
+
+        if (fnd == ht.end() || (fnd() != iter()))
+        {
+            return false;
+        }
+    }
+
+    // Are all ht elements in me?
+    for (const_iterator iter = ht.begin(); iter != ht.end(); ++iter)
+    {
+        const_iterator fnd = find(iter.key());
+
+        if (fnd == end() || (fnd() != iter()))
+        {
+            return false;
+        }
+    }
+    return true;
+}
+
+
+template<class T, class Key, class Hash>
+bool HashTable<T, Key, Hash>::operator!=(const HashTable<T, Key, Hash>& ht)
+ const
+{
+    return !(operator==(ht));
+}
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam

src/OpenFOAM/containers/HashTables/HashTable/HashTable.H

@@ -223,6 +223,15 @@ public:
         //- Assignment
         void operator=(const HashTable<T, Key, Hash>&);
 
+        //- Equality. Two hashtables are equal if all contents of first are
+        //  also in second and vice versa. So does not depend on table size or
+        //  order!
+        bool operator==(const HashTable<T, Key, Hash>&) const;
+
+        //- The opposite of the equality operation. Takes linear time.
+        bool operator!=(const HashTable<T, Key, Hash>&) const;
+
+
 
     // STL type definitions
 

src/dynamicMesh/layerAdditionRemoval/removeCellLayer.C

@@ -235,7 +235,18 @@ void Foam::layerAdditionRemoval::removeCellLayer
                     mesh.faceZones()[modifiedFaceZone].whichFace(curFaceID)
                 ];
         }
-            
+        
+        label newNei;
+        
+        if (curFaceID < mesh.nInternalFaces())
+        {
+            newNei = nei[curFaceID];
+        }
+        else
+        {
+            newNei = -1;
+        }
+
         // Modify the face
         ref.setAction
         (
@@ -244,7 +255,7 @@ void Foam::layerAdditionRemoval::removeCellLayer
                 newFace,                // modified face
                 curFaceID,              // label of face being modified
                 own[curFaceID],         // owner
-                nei[curFaceID],         // neighbour
+                newNei,                 // neighbour
                 false,                  // face flip
                 mesh.boundaryMesh().whichPatch(curFaceID),// patch for face
                 false,                  // remove from zone

src/dynamicMesh/polyTopoChange/polyTopoChange/removeCells.C

@@ -214,6 +214,18 @@ void Foam::removeCells::setRefinement
     {
         label patchI = exposedPatchIDs[i];
 
+        if (patchI < 0 || patchI >= patches.size())
+        {
+            FatalErrorIn
+            (
+                "removeCells::setRefinement(const labelList&"
+                ", const labelList&, const labelList&, polyTopoChange&)"
+            )   << "Invalid patch " << patchI
+                << " for exposed face " << exposedFaceLabels[i] << endl
+                << "Valid patches 0.." << patches.size()-1
+                << abort(FatalError);
+        }
+
         if (patches[patchI].coupled())
         {
             FatalErrorIn