Merge branch 'feature-patchProbes' into 'develop'

ENH: patchProbes: allow wildcards for patches

Also write out snapped location.

Probe location output kept as in Foundation

See merge request !11

src/sampling/probes/patchProbes.C

@@ -46,32 +46,41 @@ void Foam::patchProbes::findElements(const fvMesh& mesh)
 
     const polyBoundaryMesh& bm = mesh.boundaryMesh();
 
-    label patchI = bm.findPatchID(patchName_);
+    // All the info for nearest. Construct to miss
-
-    if (patchI == -1)
-    {
-        FatalErrorIn
-        (
-            " Foam::patchProbes::findElements(const fvMesh&)"
-        )   << " Unknown patch name "
-            << patchName_ << endl
-            << exit(FatalError);
-    }
-
-     // All the info for nearest. Construct to miss
     List<mappedPatchBase::nearInfo> nearest(this->size());
 
-    const polyPatch& pp = bm[patchI];
 
-    if (pp.size() > 0)
+    const labelList patchIDs(bm.patchSet(patchNames_).sortedToc());
+
+    label nFaces = 0;
+    forAll(patchIDs, i)
     {
-        labelList bndFaces(pp.size());
+        nFaces += bm[patchIDs[i]].size();
-        forAll(bndFaces, i)
+    }
+
+    if (nFaces > 0)
+    {
+        // Collect mesh faces and bounding box
+        labelList bndFaces(nFaces);
+        treeBoundBox overallBb(treeBoundBox::invertedBox);
+
+        nFaces = 0;
+        forAll(patchIDs, i)
         {
-            bndFaces[i] =  pp.start() + i;
+            const polyPatch& pp = bm[patchIDs[i]];
+            forAll(pp, i)
+            {
+                bndFaces[nFaces++] = pp.start()+i;
+                const face& f = pp[i];
+                forAll(f, fp)
+                {
+                    const point& pt = pp.points()[f[fp]];
+                    overallBb.min() = min(overallBb.min(), pt);
+                    overallBb.max() = max(overallBb.max(), pt);
+                }
+            }
         }
 
-        treeBoundBox overallBb(pp.points());
         Random rndGen(123456);
         overallBb = overallBb.extend(rndGen, 1e-4);
         overallBb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
@@ -106,11 +115,7 @@ void Foam::patchProbes::findElements(const fvMesh& mesh)
 
             if (!info.hit())
             {
-                info = boundaryTree.findNearest
+                info = boundaryTree.findNearest(sample, Foam::sqr(GREAT));
-                (
-                    sample,
-                    Foam::sqr(GREAT)
-                );
             }
 
             label faceI = boundaryTree.shapes().faceLabels()[info.index()];
@@ -129,20 +134,26 @@ void Foam::patchProbes::findElements(const fvMesh& mesh)
                 << " This sample will not be included "
                 << endl;
             }
-            else
+            else if (info.hit())
             {
-                const point& fc = mesh.faceCentres()[faceI];
+                // Note: do we store the face centre or the actual nearest?
+                // We interpolate using the faceI only though (no
+                // interpolation) so it does not actually matter much, just for
+                // the location written to the header.
+
+                //const point& facePt = mesh.faceCentres()[faceI];
+                const point& facePt = info.hitPoint();
 
                 mappedPatchBase::nearInfo sampleInfo;
 
                 sampleInfo.first() = pointIndexHit
                 (
                     true,
-                    fc,
+                    facePt,
                     faceI
                 );
 
-                sampleInfo.second().first() = magSqr(fc-sample);
+                sampleInfo.second().first() = magSqr(facePt-sample);
                 sampleInfo.second().second() = Pstream::myProcNo();
 
                 nearest[probeI]= sampleInfo;
@@ -155,6 +166,14 @@ void Foam::patchProbes::findElements(const fvMesh& mesh)
     Pstream::listCombineGather(nearest, mappedPatchBase::nearestEqOp());
     Pstream::listCombineScatter(nearest);
 
+
+    // Update actual probe locations
+    forAll(nearest, sampleI)
+    {
+        operator[](sampleI) = nearest[sampleI].first().rawPoint();
+    }
+
+
     if (debug)
     {
         Info<< "patchProbes::findElements" << " : " << endl;
@@ -165,26 +184,41 @@ void Foam::patchProbes::findElements(const fvMesh& mesh)
 
             Info<< "    " << sampleI << " coord:"<< operator[](sampleI)
                 << " found on processor:" << procI
-                << " in local cell/face:" << localI
+                << " in local face:" << localI
-                << " with fc:" << nearest[sampleI].first().rawPoint() << endl;
+                << " with location:" << nearest[sampleI].first().rawPoint()
+                << endl;
         }
     }
 
 
     // Extract any local faces to sample
-    elementList_.setSize(nearest.size(), -1);
+    elementList_.setSize(nearest.size());
+    elementList_ = -1;
+    faceList_.setSize(nearest.size());
+    faceList_ = -1;
 
     forAll(nearest, sampleI)
     {
         if (nearest[sampleI].second().second() == Pstream::myProcNo())
         {
             // Store the face to sample
-            elementList_[sampleI] = nearest[sampleI].first().index();
+            faceList_[sampleI] = nearest[sampleI].first().index();
         }
     }
 }
 
 
+void Foam::patchProbes::readDict(const dictionary& dict)
+{
+    if (!dict.readIfPresent("patches", patchNames_))
+    {
+        word patchName(dict.lookup("patchName"));
+        patchNames_ = wordReList(1, wordRe(patchName));
+    }
+    probes::readDict(dict);
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::patchProbes::patchProbes
@@ -192,19 +226,22 @@ Foam::patchProbes::patchProbes
     const word& name,
     const objectRegistry& obr,
     const dictionary& dict,
-    const bool loadFromFiles
+    const bool loadFromFiles,
+    const bool doFindElements
 )
 :
-    probes(name, obr, dict, loadFromFiles)
+    probes(name, obr, dict, loadFromFiles, false)
 {
-    // When constructing probes above it will have called the
+    readDict(dict);
-    // probes::findElements (since the virtual mechanism not yet operating).
-    // Not easy to workaround (apart from feeding through flag into constructor)
-    // so clear out any cells found for now.
-    elementList_.clear();
-    faceList_.clear();
 
-    read(dict);
+    if (doFindElements)
+    {
+        // Find the elements
+        findElements(mesh_);
+
+        // Open the probe streams
+        prepare();
+    }
 }
 
 
@@ -232,10 +269,16 @@ void Foam::patchProbes::write()
     }
 }
 
+
 void Foam::patchProbes::read(const dictionary& dict)
 {
-    dict.lookup("patchName") >> patchName_;
+    readDict(dict);
-    probes::read(dict);
+
+    // Find the elements
+    findElements(mesh_);
+
+    // Open the probe streams
+    prepare();
 }
 
 

src/sampling/probes/patchProbes.H

@@ -25,9 +25,45 @@ Class
     Foam::patchProbes
 
 Description
-    Set of locations to sample.at patches
+    Set of locations to sample at patches
 
     Call write() to sample and write files.
+    - find nearest location on nearest face
+    - update *this with location (so header contains 'snapped' locations
+    - use *this as the sampling location
+
+    Example of function object specification:
+    \verbatim
+    patchProbes
+    {
+        type            patchProbes;
+        functionObjectLibs ( "libsampling.so" );
+
+        // Name of the directory for probe data
+        name            patchProbes;
+
+        // Patches to sample (wildcards allowed)
+        patches         (".*inl.*");
+
+        // Write at same frequency as fields
+        outputControl   outputTime;
+        outputInterval  1;
+
+        // Fields to be probed
+        fields
+        (
+            p U
+        );
+
+        // Locations to probe. These get snapped onto the nearest point
+        // on the selected patches
+        probeLocations
+        (
+            ( -100 0 0.01 )      // at inlet
+        );
+    }
+    \endverbatim
+
 
 SourceFiles
     patchProbes.C
@@ -58,13 +94,15 @@ class patchProbes
 :
     public probes
 {
-    // Private data
+protected:
 
-        //- Patch name
+    // Protected data
-        word patchName_;
+
+        //- Patches to sample
+        wordReList patchNames_;
 
 
-    // Private Member Functions
+    // Protected Member Functions
 
         //- Sample and write a particular volume field
         template<class Type>
@@ -73,7 +111,6 @@ class patchProbes
             const GeometricField<Type, fvPatchField, volMesh>&
         );
 
-
         //- Sample and write a particular surface field
         template<class Type>
         void sampleAndWrite
@@ -81,17 +118,14 @@ class patchProbes
             const GeometricField<Type, fvsPatchField, surfaceMesh>&
         );
 
-
         //- Sample and write all the fields of the given type
         template<class Type>
         void sampleAndWrite(const fieldGroup<Type>&);
 
-
          //- Sample and write all the surface fields of the given type
         template<class Type>
         void sampleAndWriteSurfaceFields(const fieldGroup<Type>&);
 
-
         //- Sample a volume field at all locations
         template<class Type>
         tmp<Field<Type> > sample
@@ -99,7 +133,6 @@ class patchProbes
             const GeometricField<Type, fvPatchField, volMesh>&
         ) const;
 
-
         //- Sample a surface field at all locations
         template<class Type>
         tmp<Field<Type> > sample
@@ -107,11 +140,18 @@ class patchProbes
             const GeometricField<Type, fvsPatchField, surfaceMesh>&
         ) const;
 
-
         //- Sample a single field on all sample locations
         template<class Type>
         tmp<Field<Type> > sample(const word& fieldName) const;
 
+        //- Find elements containing patchProbes
+        virtual void findElements(const fvMesh&);
+
+        //- Read dictionary settings
+        void readDict(const dictionary& dict);
+
+
+private:
 
         //- Disallow default bitwise copy construct
         patchProbes(const patchProbes&);
@@ -135,7 +175,8 @@ public:
             const word& name,
             const objectRegistry&,
             const dictionary&,
-            const bool loadFromFiles = false
+            const bool loadFromFiles = false,
+            const bool findElements = true
         );
 
 
@@ -149,9 +190,6 @@ public:
 
         //- Read
         virtual void read(const dictionary&);
-
-        //- Find elements containing patchProbes
-        virtual void findElements(const fvMesh&);
 };
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/sampling/probes/patchProbesTemplates.C

@@ -208,7 +208,7 @@ Foam::patchProbes::sample
 
     forAll(*this, probeI)
     {
-        label faceI = elementList_[probeI];
+        label faceI = faceList_[probeI];
 
         if (faceI >= 0)
         {
@@ -259,7 +259,7 @@ Foam::patchProbes::sample
 
     forAll(*this, probeI)
     {
-        label faceI = elementList_[probeI];
+        label faceI = faceList_[probeI];
 
         if (faceI >= 0)
         {
@@ -274,4 +274,6 @@ Foam::patchProbes::sample
 
     return tValues;
 }
+
+
 // ************************************************************************* //

src/sampling/probes/probes.C

@@ -265,6 +265,25 @@ Foam::label Foam::probes::prepare()
 }
 
 
+void Foam::probes::readDict(const dictionary& dict)
+{
+    dict.lookup("probeLocations") >> *this;
+    dict.lookup("fields") >> fieldSelection_;
+
+    dict.readIfPresent("fixedLocations", fixedLocations_);
+    if (dict.readIfPresent("interpolationScheme", interpolationScheme_))
+    {
+        if (!fixedLocations_ && interpolationScheme_ != "cell")
+        {
+            WarningIn("void Foam::probes::read(const dictionary&)")
+                << "Only cell interpolation can be applied when "
+                << "not using fixedLocations. InterpolationScheme "
+                << "entry will be ignored";
+        }
+    }
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::probes::probes
@@ -272,7 +291,8 @@ Foam::probes::probes
     const word& name,
     const objectRegistry& obr,
     const dictionary& dict,
-    const bool loadFromFiles
+    const bool loadFromFiles,
+    const bool doFindElements
 )
 :
     pointField(0),
@@ -283,7 +303,18 @@ Foam::probes::probes
     fixedLocations_(true),
     interpolationScheme_("cell")
 {
-    read(dict);
+    // Read dictionary (but do not search for elements)
+    readDict(dict);
+
+    // Optionally find elements in constructor
+    if (doFindElements)
+    {
+        // Find the elements
+        findElements(mesh_);
+
+        // Open the probe streams
+        prepare();
+    }
 }
 
 
@@ -334,24 +365,12 @@ void Foam::probes::write()
 
 void Foam::probes::read(const dictionary& dict)
 {
-    dict.lookup("probeLocations") >> *this;
+    readDict(dict);
-    dict.lookup("fields") >> fieldSelection_;
 
-    dict.readIfPresent("fixedLocations", fixedLocations_);
+    // Find the elements
-    if (dict.readIfPresent("interpolationScheme", interpolationScheme_))
-    {
-        if (!fixedLocations_ && interpolationScheme_ != "cell")
-        {
-            WarningIn("void Foam::probes::read(const dictionary&)")
-                << "Only cell interpolation can be applied when "
-                << "not using fixedLocations.  InterpolationScheme "
-                << "entry will be ignored";
-        }
-    }
-
-    // Initialise cells to sample from supplied locations
     findElements(mesh_);
 
+    // Open the probe streams
     prepare();
 }
 
@@ -382,20 +401,28 @@ void Foam::probes::updateMesh(const mapPolyMesh& mpm)
             forAll(elementList_, i)
             {
                 label cellI = elementList_[i];
-                label newCellI = reverseMap[cellI];
+                if (cellI != -1)
-                if (newCellI == -1)
                 {
-                    // cell removed
+                    label newCellI = reverseMap[cellI];
-                }
+                    if (newCellI == -1)
-                else if (newCellI < -1)
+                    {
-                {
+                        // cell removed
-                    // cell merged
+                    }
-                    elems.append(-newCellI - 2);
+                    else if (newCellI < -1)
+                    {
+                        // cell merged
+                        elems.append(-newCellI - 2);
+                    }
+                    else
+                    {
+                        // valid new cell
+                        elems.append(newCellI);
+                    }
                 }
                 else
                 {
-                    // valid new cell
+                    // Keep -1 elements so the size stays the same
-                    elems.append(newCellI);
+                    elems.append(-1);
                 }
             }
 
@@ -410,20 +437,28 @@ void Foam::probes::updateMesh(const mapPolyMesh& mpm)
             forAll(faceList_, i)
             {
                 label faceI = faceList_[i];
-                label newFaceI = reverseMap[faceI];
+                if (faceI != -1)
-                if (newFaceI == -1)
                 {
-                    // face removed
+                    label newFaceI = reverseMap[faceI];
-                }
+                    if (newFaceI == -1)
-                else if (newFaceI < -1)
+                    {
-                {
+                        // face removed
-                    // face merged
+                    }
-                    elems.append(-newFaceI - 2);
+                    else if (newFaceI < -1)
+                    {
+                        // face merged
+                        elems.append(-newFaceI - 2);
+                    }
+                    else
+                    {
+                        // valid new face
+                        elems.append(newFaceI);
+                    }
                 }
                 else
                 {
-                    // valid new face
+                    // Keep -1 elements
-                    elems.append(newFaceI);
+                    elems.append(-1);
                 }
             }
 

src/sampling/probes/probes.H

@@ -32,6 +32,42 @@ Description
 
     Call write() to sample and write files.
 
+    Example of function object specification:
+    \verbatim
+    probes
+    {
+        type            probes;
+        functionObjectLibs ( "libsampling.so" );
+
+        // Name of the directory for probe data
+        name            probes;
+
+        // Write at same frequency as fields
+        outputControl   outputTime;
+        outputInterval  1;
+
+        // Fields to be probed
+        fields
+        (
+            p U
+        );
+
+        // Optional: do not recalculate cells if mesh moves
+        fixedLocations  false;
+
+        // Optional: interpolation scheme to use (default is cell)
+        interpolationScheme cellPoint;
+
+        probeLocations
+        (
+            ( 1e-06 0 0.01 )      // at inlet
+            (0.21 -0.20999 0.01)  // at outlet1
+            (0.21 0.20999 0.01)   // at outlet2
+            (0.21 0 0.01)         // at central block
+        );
+    }
+    \endverbatim
+
 SourceFiles
     probes.C
 
@@ -91,13 +127,13 @@ protected:
 
         //- Name of this set of probes,
         //  Also used as the name of the probes directory.
-        word name_;
+        const word name_;
 
         //- Const reference to fvMesh
         const fvMesh& mesh_;
 
         //- Load fields from files (not from objectRegistry)
-        bool loadFromFiles_;
+        const bool loadFromFiles_;
 
 
         // Read from dictonary
@@ -106,7 +142,7 @@ protected:
             wordReList fieldSelection_;
 
             //- Fixed locations, default = yes
-            //  Note: set to false for moving mesh calations where locations
+            //  Note: set to false for moving mesh calculations where locations
             //        should move with the mesh
             bool fixedLocations_;
 
@@ -141,7 +177,7 @@ protected:
             HashPtrTable<OFstream> probeFilePtrs_;
 
 
-    // Private Member Functions
+    // Protected Member Functions
 
         //- Clear old field groups
         void clearFieldGroups();
@@ -159,6 +195,8 @@ protected:
         //  returns number of fields to sample
         label prepare();
 
+        //- Read dictionary settings
+        void readDict(const dictionary& dict);
 
 private:
 
@@ -207,7 +245,8 @@ public:
             const word& name,
             const objectRegistry&,
             const dictionary&,
-            const bool loadFromFiles = false
+            const bool loadFromFiles = false,
+            const bool findElements = true
         );