ENH: support surfaceFieldValue on multiple faceZones or patches (#1874)

- additional "names" entry to specify a word/regex list of selections
  For example,
  {
      type    patch;
      name    inlets;
      names   ("inlet_[0-9].*" inlet);
  }

- if "names" exists AND contains a literal (non-regex) that can be used
  as a suitable value for "name", the "name" entry becomes optional.
  For example,
  {
      type    patch;
      names   ("inlet_[0-9].*" inlet);

      // inferred name = inlet
  }

- reduce some overhead in surfaceFieldValue

TUT: surfaceFieldValue on patches : reactingParcelFoam/verticalChannel
This commit is contained in:
Mark Olesen
2020-10-13 15:48:15 +02:00
parent f82af7cb2a
commit 6cf8151817
7 changed files with 374 additions and 167 deletions

View File

@ -32,7 +32,7 @@ License
#include "coupledPolyPatch.H"
#include "sampledSurface.H"
#include "mergePoints.H"
#include "indirectPrimitivePatch.H"
#include "uindirectPrimitivePatch.H"
#include "PatchTools.H"
#include "addToRunTimeSelectionTable.H"
@ -130,115 +130,225 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::obr() const
void Foam::functionObjects::fieldValues::surfaceFieldValue::setFaceZoneFaces()
{
const label zoneId = mesh_.faceZones().findZoneID(regionName_);
// Indices for all matches, already sorted
const labelList zoneIds
(
mesh_.faceZones().indices(selectionNames_)
);
if (zoneId < 0)
// Total number of faces selected
label numFaces = 0;
for (const label zoneId : zoneIds)
{
FatalErrorInFunction
<< type() << " " << name() << ": "
<< regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
<< " Unknown face zone name: " << regionName_
<< ". Valid face zones are: " << mesh_.faceZones().names()
<< nl << exit(FatalError);
numFaces += mesh_.faceZones()[zoneId].size();
}
const faceZone& fZone = mesh_.faceZones()[zoneId];
DynamicList<label> faceIds(fZone.size());
DynamicList<label> facePatchIds(fZone.size());
DynamicList<bool> faceFlip(fZone.size());
forAll(fZone, i)
if (zoneIds.empty())
{
const label facei = fZone[i];
FatalErrorInFunction
<< type() << ' ' << name() << ": "
<< regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
<< " No matching face zone(s): "
<< flatOutput(selectionNames_) << nl
<< " Known face zones: "
<< flatOutput(mesh_.faceZones().names()) << nl
<< exit(FatalError);
}
label faceId = -1;
label facePatchId = -1;
if (mesh_.isInternalFace(facei))
// Could also check this
#if 0
if (!returnReduce(bool(numFaces), orOp<bool>()))
{
WarningInFunction
<< type() << ' ' << name() << ": "
<< regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
<< " The faceZone specification: "
<< flatOutput(selectionNames_) << nl
<< " resulted in 0 faces" << nl
<< exit(FatalError);
}
#endif
faceId_.resize(numFaces);
facePatchId_.resize(numFaces);
faceFlip_.resize(numFaces);
numFaces = 0;
for (const label zoneId : zoneIds)
{
const faceZone& fZone = mesh_.faceZones()[zoneId];
forAll(fZone, i)
{
faceId = facei;
facePatchId = -1;
}
else
{
facePatchId = mesh_.boundaryMesh().whichPatch(facei);
const polyPatch& pp = mesh_.boundaryMesh()[facePatchId];
if (isA<coupledPolyPatch>(pp))
const label meshFacei = fZone[i];
const bool isFlip = fZone.flipMap()[i];
// Internal faces
label faceId = meshFacei;
label facePatchId = -1;
// Boundary faces
if (!mesh_.isInternalFace(meshFacei))
{
if (refCast<const coupledPolyPatch>(pp).owner())
facePatchId = mesh_.boundaryMesh().whichPatch(meshFacei);
const polyPatch& pp = mesh_.boundaryMesh()[facePatchId];
if (isA<coupledPolyPatch>(pp))
{
faceId = pp.whichFace(facei);
if (refCast<const coupledPolyPatch>(pp).owner())
{
faceId = pp.whichFace(meshFacei);
}
else
{
faceId = -1;
}
}
else if (!isA<emptyPolyPatch>(pp))
{
faceId = meshFacei - pp.start();
}
else
{
faceId = -1;
facePatchId = -1;
}
}
else if (!isA<emptyPolyPatch>(pp))
{
faceId = facei - pp.start();
}
else
{
faceId = -1;
facePatchId = -1;
}
}
if (faceId >= 0)
{
faceIds.append(faceId);
facePatchIds.append(facePatchId);
faceFlip.append(fZone.flipMap()[i] ? true : false);
if (faceId >= 0)
{
faceId_[numFaces] = faceId;
facePatchId_[numFaces] = facePatchId;
faceFlip_[numFaces] = isFlip;
++numFaces;
}
}
}
faceId_.transfer(faceIds);
facePatchId_.transfer(facePatchIds);
faceFlip_.transfer(faceFlip);
// Shrink to size used
faceId_.resize(numFaces);
facePatchId_.resize(numFaces);
faceFlip_.resize(numFaces);
nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
if (debug)
{
Pout<< "Original face zone size = " << fZone.size()
<< ", new size = " << faceId_.size() << endl;
}
}
void Foam::functionObjects::fieldValues::surfaceFieldValue::setPatchFaces()
{
const label patchid = mesh_.boundaryMesh().findPatchID(regionName_);
// Patch indices for all matches
labelList patchIds;
if (patchid < 0)
// Total number of faces selected
label numFaces = 0;
labelList selected
(
mesh_.boundaryMesh().patchSet
(
selectionNames_,
false // warnNotFound - we do that ourselves
).sortedToc()
);
DynamicList<label> bad;
for (const label patchi : selected)
{
const polyPatch& pp = mesh_.boundaryMesh()[patchi];
if (isA<emptyPolyPatch>(pp))
{
bad.append(patchi);
}
else
{
numFaces += pp.size();
}
}
if (bad.size())
{
label nGood = (selected.size() - bad.size());
auto& os = (nGood > 0 ? WarningInFunction : FatalErrorInFunction);
os << "Cannot sample an empty patch" << nl;
for (const label patchi : bad)
{
os << " "
<< mesh_.boundaryMesh()[patchi].name() << nl;
}
if (nGood)
{
os << "No non-empty patches selected" << endl
<< exit(FatalError);
}
else
{
os << "Selected " << nGood << " non-empty patches" << nl;
}
patchIds.resize(nGood);
nGood = 0;
for (const label patchi : selected)
{
if (!bad.found(patchi))
{
patchIds[nGood] = patchi;
++nGood;
}
}
}
else
{
patchIds = std::move(selected);
}
if (patchIds.empty())
{
FatalErrorInFunction
<< type() << " " << name() << ": "
<< type() << ' ' << name() << ": "
<< regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
<< " Unknown patch name: " << regionName_
<< ". Valid patch names are: "
<< " No matching patch name(s): "
<< flatOutput(selectionNames_) << nl
<< " Known patch names:" << nl
<< mesh_.boundaryMesh().names() << nl
<< exit(FatalError);
}
const polyPatch& pp = mesh_.boundaryMesh()[patchid];
label nFaces = pp.size();
if (isA<emptyPolyPatch>(pp))
// Could also check this
#if 0
if (!returnReduce(bool(numFaces), orOp<bool>()))
{
nFaces = 0;
WarningInFunction
<< type() << ' ' << name() << ": "
<< regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
<< " The patch specification: "
<< flatOutput(selectionNames_) << nl
<< " resulted in 0 faces" << nl
<< exit(FatalError);
}
#endif
faceId_.setSize(nFaces);
facePatchId_.setSize(nFaces);
faceFlip_.setSize(nFaces);
faceId_.resize(numFaces);
facePatchId_.resize(numFaces);
faceFlip_.resize(numFaces);
nFaces_ = returnReduce(faceId_.size(), sumOp<label>());
forAll(faceId_, facei)
numFaces = 0;
for (const label patchi : patchIds)
{
faceId_[facei] = facei;
facePatchId_[facei] = patchid;
faceFlip_[facei] = false;
const polyPatch& pp = mesh_.boundaryMesh()[patchi];
const label len = pp.size();
SubList<label>(faceId_, len, numFaces) = identity(len);
SubList<label>(facePatchId_, len, numFaces) = patchi;
SubList<bool>(faceFlip_, len, numFaces) = false;
numFaces += len;
}
}
@ -252,24 +362,25 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry
List<faceList> allFaces(Pstream::nProcs());
List<pointField> allPoints(Pstream::nProcs());
labelList globalFacesIs(faceId_);
forAll(globalFacesIs, i)
{
if (facePatchId_[i] != -1)
IndirectList<face> selectedFaces(mesh_.faces(), labelList(faceId_));
labelList& meshFaceIds = selectedFaces.addressing();
forAll(meshFaceIds, i)
{
const label patchi = facePatchId_[i];
globalFacesIs[i] += mesh_.boundaryMesh()[patchi].start();
if (patchi != -1)
{
meshFaceIds[i] += mesh_.boundaryMesh()[patchi].start();
}
}
}
// Add local faces and points to the all* lists
indirectPrimitivePatch pp
(
IndirectList<face>(mesh_.faces(), globalFacesIs),
mesh_.points()
);
allFaces[Pstream::myProcNo()] = pp.localFaces();
allPoints[Pstream::myProcNo()] = pp.localPoints();
// Add local faces and points to the all* lists
uindirectPrimitivePatch pp(selectedFaces, mesh_.points());
allFaces[Pstream::myProcNo()] = pp.localFaces();
allPoints[Pstream::myProcNo()] = pp.localPoints();
}
Pstream::gatherList(allFaces);
Pstream::gatherList(allPoints);
@ -283,53 +394,41 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::combineMeshGeometry
nPoints += allPoints[proci].size();
}
faces.setSize(nFaces);
points.setSize(nPoints);
faces.resize(nFaces);
points.resize(nPoints);
nFaces = 0;
nPoints = 0;
// My own data first
// My data first
{
const faceList& fcs = allFaces[Pstream::myProcNo()];
for (const face& f : fcs)
for (const face& f : allFaces[Pstream::myProcNo()])
{
face& newF = faces[nFaces++];
newF.setSize(f.size());
forAll(f, fp)
{
newF[fp] = f[fp] + nPoints;
}
faces[nFaces++] = offsetOp<face>()(f, nPoints);
}
const pointField& pts = allPoints[Pstream::myProcNo()];
for (const point& pt : pts)
for (const point& p : allPoints[Pstream::myProcNo()])
{
points[nPoints++] = pt;
points[nPoints++] = p;
}
}
// Other proc data follows
forAll(allFaces, proci)
{
if (proci != Pstream::myProcNo())
if (proci == Pstream::myProcNo())
{
const faceList& fcs = allFaces[proci];
for (const face& f : fcs)
{
face& newF = faces[nFaces++];
newF.setSize(f.size());
forAll(f, fp)
{
newF[fp] = f[fp] + nPoints;
}
}
continue;
}
const pointField& pts = allPoints[proci];
for (const point& pt : pts)
{
points[nPoints++] = pt;
}
for (const face& f : allFaces[proci])
{
faces[nFaces++] = offsetOp<face>()(f, nPoints);
}
for (const point& p : allPoints[proci])
{
points[nPoints++] = p;
}
}
@ -383,11 +482,7 @@ combineSurfaceGeometry
PatchTools::gatherAndMerge
(
mergeDim,
primitivePatch
(
SubList<face>(s.faces(), s.faces().size()),
s.points()
),
primitivePatch(SubList<face>(s.faces()), s.points()),
points,
faces,
pointsMap
@ -413,11 +508,7 @@ combineSurfaceGeometry
PatchTools::gatherAndMerge
(
mergeDim,
primitivePatch
(
SubList<face>(s.faces(), s.faces().size()),
s.points()
),
primitivePatch(SubList<face>(s.faces()), s.points()),
points,
faces,
pointsMap
@ -522,7 +613,7 @@ bool Foam::functionObjects::fieldValues::surfaceFieldValue::update()
if (nFaces_ == 0)
{
FatalErrorInFunction
<< type() << " " << name() << ": "
<< type() << ' ' << name() << ": "
<< regionTypeNames_[regionType_] << '(' << regionName_ << "):" << nl
<< " Region has no faces" << exit(FatalError);
}
@ -548,7 +639,7 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::writeFileHeader
if (canWriteHeader() && (operation_ != opNone))
{
writeCommented(os, "Region type : ");
os << regionTypeNames_[regionType_] << " " << regionName_ << endl;
os << regionTypeNames_[regionType_] << ' ' << regionName_ << nl;
writeHeaderValue(os, "Faces", nFaces_);
writeHeaderValue(os, "Area", totalArea_);
@ -570,7 +661,7 @@ void Foam::functionObjects::fieldValues::surfaceFieldValue::writeFileHeader
for (const word& fieldName : fields_)
{
os << tab << operationTypeNames_[operation_]
<< "(" << fieldName << ")";
<< '(' << fieldName << ')';
}
os << endl;
@ -821,9 +912,10 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::surfaceFieldValue
true // Failsafe behaviour
)
),
weightFieldName_("none"),
needsUpdate_(true),
writeArea_(false),
selectionNames_(),
weightFieldName_("none"),
totalArea_(0),
nFaces_(0),
faceId_(),
@ -854,9 +946,10 @@ Foam::functionObjects::fieldValues::surfaceFieldValue::surfaceFieldValue
true // Failsafe behaviour
)
),
weightFieldName_("none"),
needsUpdate_(true),
writeArea_(false),
selectionNames_(),
weightFieldName_("none"),
totalArea_(0),
nFaces_(0),
faceId_(),
@ -876,18 +969,54 @@ bool Foam::functionObjects::fieldValues::surfaceFieldValue::read
{
fieldValue::read(dict);
weightFieldName_ = "none";
needsUpdate_ = true;
writeArea_ = dict.getOrDefault("writeArea", false);
weightFieldName_ = "none";
totalArea_ = 0;
nFaces_ = 0;
faceId_.clear();
facePatchId_.clear();
faceFlip_.clear();
sampledPtr_.clear();
surfaceWriterPtr_.clear();
sampledPtr_.reset(nullptr);
surfaceWriterPtr_.reset(nullptr);
// Can have "name" (word) and/or "names" (wordRes)
//
// If "names" exists AND contains a literal (non-regex) that can be used
// as a suitable value for "name", the "name" entry becomes optional.
regionName_.clear();
selectionNames_.clear();
{
dict.readIfPresent("names", selectionNames_);
for (const auto& item : selectionNames_)
{
if (item.isLiteral())
{
regionName_ = item;
break;
}
}
// Mandatory if we didn't pick up a value from selectionNames_
dict.readEntry
(
"name",
regionName_,
keyType::LITERAL,
regionName_.empty()
);
// Ensure there is always content for selectionNames_
if (selectionNames_.empty())
{
selectionNames_.resize(1);
selectionNames_.first() = regionName_;
}
}
dict.readEntry("name", regionName_);
// Create sampled surface, but leave 'expired' (ie, no update) since it
// may depend on fields or data that do not yet exist
@ -901,7 +1030,7 @@ bool Foam::functionObjects::fieldValues::surfaceFieldValue::read
);
}
Info<< type() << " " << name() << ":" << nl
Info<< type() << ' ' << name() << ':' << nl
<< " operation = ";
if (postOperation_ != postOpNone)

View File

@ -31,17 +31,17 @@ Group
grpFieldFunctionObjects
Description
Provides a 'face regionType' variant of the \c fieldValues function object.
A \c face regionType variant of the \c fieldValues function object.
Given a list of user-specified fields and a selection of mesh (or general
surface) faces, a number of operations can be performed, such as sums,
averages and integrations.
For example, to calculate the volumetric or mass flux across a patch,
apply the 'sum' operator to the flux field (typically \c phi).
apply the \c sum operator to the flux field (typically \c phi).
Usage
Minimal example by using \c system/controlDict.functions:
A minimal example:
\verbatim
surfaceFieldValuePatch1
{
@ -56,6 +56,8 @@ Usage
name <patch>;
// Optional entries (runtime modifiable)
names (<patch-name> <patch-regex>);
postOperation none;
weightField alpha1;
scaleFactor 1.0;
@ -79,6 +81,8 @@ Usage
name <faceZone>;
// Optional entries (runtime modifiable)
names (<zone-name> <zone-regex>);
postOperation none;
weightField alpha1;
scaleFactor 1.0;
@ -92,12 +96,13 @@ Usage
where the entries mean:
\table
Property | Description | Type | Req'd | Dflt
Property | Description | Type | Reqd | Dflt
type | Type name: surfaceFieldValue | word | yes | -
libs | Library name: fieldFunctionObjects | word | yes | -
fields | Names of operand fields | wordList | yes | -
libs | Libraries: fieldFunctionObjects | wordList | yes | -
regionType | Face regionType: see below | word | yes | -
name | Name for regionType | word | yes | -
fields | Names of operand fields | wordList | yes | -
name | Name of the regionType | word | yes | -
names | Extended selection | word/regex list | no | -
operation | Operation type: see below | word | yes | -
postOperation | Post-operation type: see below | word | no | none
weightField | Name of field to apply weighting | word | no | none
@ -112,9 +117,9 @@ Usage
Options for the \c regionType entry:
\plaintable
faceZone | The \b name entry to specify the faceZone
patch | The \b name entry to specify the patch
functionObjectSurface | The \b name entry to specify a polySurface
faceZone | The \b name entry specifies a faceZone. Supports \b names
patch | The \b name entry specifies a patch. Supports \b names
functionObjectSurface | The \b name entry specifies a polySurface
sampledSurface | A \b sampledSurfaceDict sub-dictionary and \b name
\endplaintable
@ -156,12 +161,16 @@ Usage
Usage by the \c postProcess utility is not available.
Note
- Some types (eg, patch or faceZone) support the selection of multiple
entries, which can be specified as list of words or regular expressions
by \b names entry.<br>
If the \b names enty exists \em and contains a literal that can be used
as a suitable value for \b name, the \b name entry becomes optional
instead of being mandatory.
- The values reported by the \c areaNormalAverage and \c areaNormalIntegrate
operations are written as the first component of a field with the same
rank as the input field.
- Faces on empty patches get ignored.
- If the field is a volField the \c faceZone
can only consist of boundary faces.
- Faces on empty patches are ignored.
- Using \c functionObjectSurface:
- The keyword %subRegion should not be used to select surfaces.
Instead specify the regionType 'functionObjectSurface' and provide
@ -242,8 +251,8 @@ public:
//- Region type enumeration
enum regionTypes
{
stFaceZone = 0x01, //!< Calculate on a faceZone
stPatch = 0x02, //!< Calculate on a patch
stFaceZone = 0x01, //!< Calculate with faceZone(s)
stPatch = 0x02, //!< Calculate with patch(es)
stObject = 0x11, //!< Calculate with function object surface
stSampled = 0x12 //!< Sample onto surface and calculate
};
@ -385,15 +394,18 @@ protected:
//- Optional post-evaluation operation
postOperationType postOperation_;
//- Weight field name - optional
word weightFieldName_;
//- Track if the surface needs an update
bool needsUpdate_;
//- Optionally write the area of the surfaceFieldValue
bool writeArea_;
//- Extended selections
wordRes selectionNames_;
//- Weight field name - optional
word weightFieldName_;
//- Total area of the surfaceFieldValue
scalar totalArea_;
@ -401,17 +413,20 @@ protected:
label nFaces_;
// If operating on mesh faces (faceZone, patch)
// If operating on mesh faces (faceZone, patch)
//- Local list of face IDs
labelList faceId_;
//- Local list of face IDs
labelList faceId_;
//- Local list of patch ID per face
labelList facePatchId_;
//- Local list of patch ID per face
labelList facePatchId_;
//- List representing the face flip map
// (false: use as-is, true: negate)
boolList faceFlip_;
//- List representing the face flip map
// (false: use as-is, true: negate)
boolList faceFlip_;
// Demand-driven
//- The sampledSurface (when operating on sampledSurface)
autoPtr<sampledSurface> sampledPtr_;

View File

@ -3,15 +3,12 @@ cd "${0%/*}" || exit # Run from this directory
. ${WM_PROJECT_DIR:?}/bin/tools/RunFunctions # Tutorial run functions
#------------------------------------------------------------------------------
# create mesh
runApplication blockMesh
restore0Dir
# initialise with potentialFoam solution
runApplication blockMesh
runApplication potentialFoam
# run the solver
runApplication $(getApplication)
#------------------------------------------------------------------------------

View File

@ -0,0 +1,16 @@
#!/bin/sh
cd "${0%/*}" || exit # Run from this directory
. ${WM_PROJECT_DIR:?}/bin/tools/RunFunctions # Tutorial run functions
#------------------------------------------------------------------------------
restore0Dir
runApplication blockMesh
runApplication potentialFoam
runApplication decomposePar
runParallel $(getApplication)
#------------------------------------------------------------------------------

View File

@ -53,7 +53,27 @@ maxDeltaT 1e-03;
functions
{
surfaceFieldValue1
avgInlets
{
type surfaceFieldValue;
libs (fieldFunctionObjects);
enabled yes;
writeControl writeTime;
log yes;
writeFields off;
surfaceFormat vtk;
regionType patch;
// name inlet;
names (inlet "inlet.*");
operation weightedAverage;
weightField phi;
fields
(
T
);
}
avgOutlets
{
type surfaceFieldValue;
libs (fieldFunctionObjects);

View File

@ -0,0 +1,30 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2006 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object decomposeParDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
numberOfSubdomains 4;
method scotch;
coeffs
{
n (2 2 1);
}
distributed no;
roots ( );
// ************************************************************************* //

View File

@ -48,7 +48,7 @@ runTimeModifiable yes;
functions
{
surfaceFieldValue1
avgOutlets
{
type surfaceFieldValue;
libs (fieldFunctionObjects);