ENH: cyclicPolyPatch: have local 'matchTolerance' entry in dictionary

applications/utilities/mesh/manipulation/createPatch/createPatch.C

@@ -526,13 +526,6 @@ int main(int argc, char *argv[])
     // Whether to synchronise points
     const Switch pointSync(dict.lookup("pointSync"));
 
-
-    // Set the matching tolerance so we can read illegal meshes
-    scalar tol = readScalar(dict.lookup("matchTolerance"));
-    Info<< "Using relative tolerance " << tol
-        << " to match up faces and points" << nl << endl;
-    coupledPolyPatch::matchTol = tol;
-
 #   include "createNamedPolyMesh.H"
 
     const word oldInstance = mesh.pointsInstance();

applications/utilities/mesh/manipulation/createPatch/createPatchDict

@@ -34,16 +34,12 @@ FoamFile
 // This will usually fail upon loading:
 //  "face 0 area does not match neighbour 2 by 0.0100005%"
 //  " -- possible face ordering problem."
-// - change patch type from 'cyclic' to 'patch' in the polyMesh/boundary file.
+// - in polyMesh/boundary file:
-// - loosen match tolerance to get case to load
+//      - loosen matchTolerance of all cyclics to get case to load
+//      - or change patch type from 'cyclic' to 'patch'
 // - regenerate cyclic as above
 
 
-// Tolerance used in matching faces. Absolute tolerance is span of
-// face times this factor. To load incorrectly matches meshes set this
-// to a higher value.
-matchTolerance 1E-3;
-
 // Do a synchronisation of coupled points after creation of any patches.
 // Note: this does not work with points that are on multiple coupled patches
 //       with transformations.
@@ -67,6 +63,12 @@ patches
             transform rotational;
             rotationAxis (1 0 0);
             rotationCentre (0 0 0);
+            // transform translational;
+            // separationVector (1 0 0);
+
+            // Optional non-default tolerance to be able to define cyclics
+            // on bad meshes
+            //matchTolerance 1E-2;
         }
 
         // How to construct: either from 'patches' or 'set'

src/OpenFOAM/meshes/polyMesh/polyPatches/basic/coupled/coupledPolyPatch.C

@@ -34,7 +34,7 @@ namespace Foam
 {
     defineTypeNameAndDebug(coupledPolyPatch, 0);
 
-    scalar coupledPolyPatch::matchTol = 1E-3;
+    const scalar coupledPolyPatch::defaultMatchTol_ = 1E-4;
 
     template<>
     const char* NamedEnum<coupledPolyPatch::transformType, 4>::names[] =
@@ -145,6 +145,7 @@ Foam::pointField Foam::coupledPolyPatch::getAnchorPoints
 
 Foam::scalarField Foam::coupledPolyPatch::calcFaceTol
 (
+    const scalar matchTol,
     const UList<face>& faces,
     const pointField& points,
     const pointField& faceCentres
@@ -401,7 +402,8 @@ Foam::coupledPolyPatch::coupledPolyPatch
     const polyBoundaryMesh& bm
 )
 :
-    polyPatch(name, size, start, index, bm)
+    polyPatch(name, size, start, index, bm),
+    matchTolerance_(defaultMatchTol_)
 {}
 
 
@@ -413,7 +415,8 @@ Foam::coupledPolyPatch::coupledPolyPatch
     const polyBoundaryMesh& bm
 )
 :
-    polyPatch(name, dict, index, bm)
+    polyPatch(name, dict, index, bm),
+    matchTolerance_(dict.lookupOrDefault("matchTolerance", defaultMatchTol_))
 {}
 
 
@@ -423,7 +426,8 @@ Foam::coupledPolyPatch::coupledPolyPatch
     const polyBoundaryMesh& bm
 )
 :
-    polyPatch(pp, bm)
+    polyPatch(pp, bm),
+    matchTolerance_(pp.matchTolerance_)
 {}
 
 
@@ -436,7 +440,8 @@ Foam::coupledPolyPatch::coupledPolyPatch
     const label newStart
 )
 :
-    polyPatch(pp, bm, index, newSize, newStart)
+    polyPatch(pp, bm, index, newSize, newStart),
+    matchTolerance_(pp.matchTolerance_)
 {}
 
 
@@ -449,7 +454,8 @@ Foam::coupledPolyPatch::coupledPolyPatch
     const label newStart
 )
 :
-    polyPatch(pp, bm, index, mapAddressing, newStart)
+    polyPatch(pp, bm, index, mapAddressing, newStart),
+    matchTolerance_(pp.matchTolerance_)
 {}
 
 
@@ -459,4 +465,17 @@ Foam::coupledPolyPatch::~coupledPolyPatch()
 {}
 
 
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+void Foam::coupledPolyPatch::write(Ostream& os) const
+{
+    polyPatch::write(os);
+    //if (matchTolerance_ != defaultMatchTol_)
+    {
+        os.writeKeyword("matchTolerance") << matchTolerance_
+            << token::END_STATEMENT << nl;
+    }
+}
+
+
 // ************************************************************************* //

src/OpenFOAM/meshes/polyMesh/polyPatches/basic/coupled/coupledPolyPatch.H

@@ -69,6 +69,12 @@ private:
 
     // Private data
 
+        //- default matching tolerance
+        static const scalar defaultMatchTol_;
+
+        //- local matching tolerance
+        const scalar matchTolerance_;
+
         //- offset (distance) vector from one side of the couple to the other
         mutable vectorField separation_;
 
@@ -81,14 +87,6 @@ private:
         //- Are faces collocated. Either size 0,1 or length of patch.
         mutable boolList collocated_;
 
-public:
-
-    // Static data members
-
-        //- Relative tolerance (for geometric matching).
-        static scalar matchTol;
-
-
 protected:
 
     // Protected Member Functions
@@ -105,7 +103,7 @@ protected:
             const vectorField& nf,
             const vectorField& nr,
             const scalarField& smallDist,
-            const scalar absTol = matchTol,
+            const scalar absTol,
             const transformType = UNKNOWN
         ) const;
 
@@ -162,6 +160,7 @@ protected:
         //  from face centre to any of the face vertices.
         static scalarField calcFaceTol
         (
+            const scalar matchTol,
             const UList<face>& faces,
             const pointField& points,
             const pointField& faceCentres
@@ -295,6 +294,11 @@ public:
                 return collocated_;
             }
 
+            scalar matchTolerance() const
+            {
+                return matchTolerance_;
+            }
+
 
         //- Calculate the patch geometry
         virtual void calcGeometry
@@ -328,6 +332,9 @@ public:
             labelList& faceMap,
             labelList& rotation
         ) const = 0;
+
+        //- Write the polyPatch data as a dictionary
+        virtual void write(Ostream&) const;
 };
 
 

src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/cyclic/cyclicPolyPatch.C

@@ -212,7 +212,7 @@ void Foam::cyclicPolyPatch::calcTransforms
                     maxAreaFacei = facei;
                 }
 
-                if (areaDiff > coupledPolyPatch::matchTol)
+                if (areaDiff > matchTolerance())
                 {
                     FatalErrorIn
                     (
@@ -224,13 +224,17 @@ void Foam::cyclicPolyPatch::calcTransforms
                         << "patch:" << name()
                         << " my area:" << magSf
                         << " neighbour area:" << nbrMagSf
-                        << " matching tolerance:" << coupledPolyPatch::matchTol
+                        << " matching tolerance:" << matchTolerance()
                          << endl
                         << "Mesh face:" << start()+facei
                         << " fc:" << half0Ctrs[facei]
                         << endl
                         << "Neighbour fc:" << half1Ctrs[facei]
                         << endl
+                        << "If you are certain your matching is correct"
+                        << " you can increase the 'matchTolerance' setting"
+                        << " in the patch dictionary in the boundary file."
+                        << endl
                         << "Rerun with cyclic debug flag set"
                         << " for more information." << exit(FatalError);
                 }
@@ -302,6 +306,7 @@ void Foam::cyclicPolyPatch::calcTransforms
             (
                 calcFaceTol
                 (
+                    matchTolerance(),
                     half0,
                     half0.points(),
                     static_cast<const pointField&>(half0Ctrs)
@@ -315,7 +320,7 @@ void Foam::cyclicPolyPatch::calcTransforms
                 half0Normals,
                 half1Normals,
                 half0Tols,
-                matchTol,
+                matchTolerance(),
                 transform_
             );
 
@@ -506,7 +511,7 @@ void Foam::cyclicPolyPatch::getCentresAndAnchors
             vector n1 = pp1[max1I].normal(pp1.points());
             n1 /= mag(n1) + VSMALL;
 
-            if (mag(n0 & n1) < 1-coupledPolyPatch::matchTol)
+            if (mag(n0 & n1) < 1-matchTolerance())
             {
                 if (debug)
                 {
@@ -557,7 +562,7 @@ void Foam::cyclicPolyPatch::getCentresAndAnchors
 
 
     // Calculate typical distance per face
-    tols = calcFaceTol(pp1, pp1.points(), half1Ctrs);
+    tols = calcFaceTol(matchTolerance(), pp1, pp1.points(), half1Ctrs);
 }
 
 
@@ -1445,7 +1450,7 @@ bool Foam::cyclicPolyPatch::order
 
 void Foam::cyclicPolyPatch::write(Ostream& os) const
 {
-    polyPatch::write(os);
+    coupledPolyPatch::write(os);
     os.writeKeyword("neighbourPatch") << neighbPatchName_
         << token::END_STATEMENT << nl;
     switch (transform_)

src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/oldCyclic/oldCyclicPolyPatch.C

@@ -386,7 +386,7 @@ void Foam::oldCyclicPolyPatch::getCentresAndAnchors
             vector n1 = half1Faces[max1I].normal(pp.points());
             n1 /= mag(n1) + VSMALL;
 
-            if (mag(n0 & n1) < 1-coupledPolyPatch::matchTol)
+            if (mag(n0 & n1) < 1-matchTolerance())
             {
                 if (debug)
                 {
@@ -444,7 +444,7 @@ void Foam::oldCyclicPolyPatch::getCentresAndAnchors
 
 
     // Calculate typical distance per face
-    tols = calcFaceTol(half1Faces, pp.points(), half1Ctrs);
+    tols = calcFaceTol(matchTolerance(), half1Faces, pp.points(), half1Ctrs);
 }
 
 

src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processor/processorPolyPatch.C

@@ -194,7 +194,7 @@ void Foam::processorPolyPatch::calcGeometry(PstreamBuffers& pBufs)
                 faceNormals[facei] = point(1, 0, 0);
                 nbrFaceNormals[facei] = faceNormals[facei];
             }
-            else if (mag(magSf - nbrMagSf)/avSf > coupledPolyPatch::matchTol)
+            else if (mag(magSf - nbrMagSf)/avSf > matchTolerance())
             {
                 fileName nm
                 (
@@ -214,12 +214,16 @@ void Foam::processorPolyPatch::calcGeometry(PstreamBuffers& pBufs)
                     << "patch:" << name()
                     << " my area:" << magSf
                     << " neighbour area:" << nbrMagSf
-                    << " matching tolerance:" << coupledPolyPatch::matchTol
+                    << " matching tolerance:" << matchTolerance()
                     << endl
                     << "Mesh face:" << start()+facei
                     << " vertices:"
                     << UIndirectList<point>(points(), operator[](facei))()
                     << endl
+                    << "If you are certain your matching is correct"
+                    << " you can increase the 'matchTolerance' setting"
+                    << " in the patch dictionary in the boundary file."
+                    << endl
                     << "Rerun with processor debug flag set for"
                     << " more information." << exit(FatalError);
             }
@@ -236,7 +240,8 @@ void Foam::processorPolyPatch::calcGeometry(PstreamBuffers& pBufs)
             neighbFaceCentres_,
             faceNormals,
             nbrFaceNormals,
-            calcFaceTol(*this, points(), faceCentres())
+            calcFaceTol(matchTolerance(), *this, points(), faceCentres()),
+            matchTolerance()
         );
     }
 }
@@ -538,7 +543,10 @@ bool Foam::processorPolyPatch::order
         }
 
         // Calculate typical distance from face centre
-        scalarField tols(calcFaceTol(pp, pp.points(), pp.faceCentres()));
+        scalarField tols
+        (
+            calcFaceTol(matchTolerance(), pp, pp.points(), pp.faceCentres())
+        );
 
         if (debug || masterCtrs.size() != pp.size())
         {
@@ -697,7 +705,7 @@ bool Foam::processorPolyPatch::order
 
 void Foam::processorPolyPatch::write(Ostream& os) const
 {
-    polyPatch::write(os);
+    coupledPolyPatch::write(os);
     os.writeKeyword("myProcNo") << myProcNo_
         << token::END_STATEMENT << nl;
     os.writeKeyword("neighbProcNo") << neighbProcNo_

src/OpenFOAM/primitives/globalIndexAndTransform/globalIndexAndTransform.C

@@ -125,6 +125,7 @@ void Foam::globalIndexAndTransform::determineTransforms()
     const polyBoundaryMesh& patches = mesh_.boundaryMesh();
 
     transforms_ = List<vectorTensorTransform>(6);
+    scalarField maxTol(6);
 
     label nextTrans = 0;
 
@@ -148,8 +149,6 @@ void Foam::globalIndexAndTransform::determineTransforms()
 
                     if (mag(sepVec) > SMALL)
                     {
-                        scalar tol = coupledPolyPatch::matchTol;
-
                         vectorTensorTransform transform(sepVec);
 
                         if
@@ -159,12 +158,13 @@ void Foam::globalIndexAndTransform::determineTransforms()
                                 transforms_,
                                 dummyMatch,
                                 transform,
-                                tol,
+                                cpp.matchTolerance(),
                                 false
                             ) == 0
                         )
                         {
-                            transforms_[nextTrans++] = transform;
+                            transforms_[nextTrans] = transform;
+                            maxTol[nextTrans++] = cpp.matchTolerance();
                         }
 
                         if (nextTrans > 6)
@@ -191,8 +191,6 @@ void Foam::globalIndexAndTransform::determineTransforms()
 
                     if (mag(transT - I) > SMALL)
                     {
-                        scalar tol = coupledPolyPatch::matchTol;
-
                         vectorTensorTransform transform(transT);
 
                         if
@@ -202,12 +200,13 @@ void Foam::globalIndexAndTransform::determineTransforms()
                                 transforms_,
                                 dummyMatch,
                                 transform,
-                                tol,
+                                cpp.matchTolerance(),
                                 false
                             ) == 0
                         )
                         {
-                            transforms_[nextTrans++] = transform;
+                            transforms_[nextTrans] = transform;
+                            maxTol[nextTrans++] = cpp.matchTolerance();
                         }
 
                         if (nextTrans > 6)
@@ -227,12 +226,18 @@ void Foam::globalIndexAndTransform::determineTransforms()
         }
     }
 
+
+    // Collect transforms on master
+
     List<List<vectorTensorTransform> > allTransforms(Pstream::nProcs());
-
     allTransforms[Pstream::myProcNo()] = transforms_;
-
     Pstream::gatherList(allTransforms);
 
+    // Collect matching tolerance on master
+    List<scalarField> allTols(Pstream::nProcs());
+    allTols[Pstream::myProcNo()] = maxTol;
+    Pstream::gatherList(allTols);
+
     if (Pstream::master())
     {
         transforms_ = List<vectorTensorTransform>(3);
@@ -250,8 +255,6 @@ void Foam::globalIndexAndTransform::determineTransforms()
 
                 if (mag(transform.t()) > SMALL || transform.hasR())
                 {
-                    scalar tol = coupledPolyPatch::matchTol;
-
                     if
                     (
                         matchTransform
@@ -259,7 +262,7 @@ void Foam::globalIndexAndTransform::determineTransforms()
                             transforms_,
                             dummyMatch,
                             transform,
-                            tol,
+                            allTols[procI][pSVI],
                             true
                         ) ==  0
                     )
@@ -378,8 +381,6 @@ void Foam::globalIndexAndTransform::determinePatchTransformSign()
 
                     if (mag(sepVec) > SMALL)
                     {
-                        scalar tol = coupledPolyPatch::matchTol;
-
                         vectorTensorTransform t(sepVec);
 
                         label sign = matchTransform
@@ -387,7 +388,7 @@ void Foam::globalIndexAndTransform::determinePatchTransformSign()
                             transforms_,
                             matchTransI,
                             t,
-                            tol,
+                            cpp.matchTolerance(),
                             true
                         );
 
@@ -424,8 +425,6 @@ void Foam::globalIndexAndTransform::determinePatchTransformSign()
 
                     if (mag(transT - I) > SMALL)
                     {
-                        scalar tol = coupledPolyPatch::matchTol;
-
                         vectorTensorTransform t(transT);
 
                         label sign = matchTransform
@@ -433,7 +432,7 @@ void Foam::globalIndexAndTransform::determinePatchTransformSign()
                             transforms_,
                             matchTransI,
                             t,
-                            tol,
+                            cpp.matchTolerance(),
                             true
                         );
 

src/finiteVolume/fvMesh/fvPatches/constraint/cyclic/cyclicFvPatch.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2004-2011 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -52,17 +52,6 @@ void Foam::cyclicFvPatch::makeWeights(scalarField& w) const
 
     forAll(magFa, facei)
     {
-        scalar avFa = (magFa[facei] + nbrMagFa[facei])/2.0;
-
-        if (mag(magFa[facei] - nbrMagFa[facei])/avFa > 1e-4)
-        {
-            FatalErrorIn("cyclicFvPatch::makeWeights(scalarField&) const")
-                << "face " << facei << " areas do not match by "
-                << 100*mag(magFa[facei] - nbrMagFa[facei])/avFa
-                << "% -- possible face ordering problem"
-                << abort(FatalError);
-        }
-
         scalar di = deltas[facei];
         scalar dni = nbrDeltas[facei];