Go to the documentation of this file.
48 void Foam::thresholdCellFaces::calculate
51 const scalar lowerThreshold,
52 const scalar upperThreshold,
53 const bool triangulate
93 DynamicList<label> surfCells(
surfFaces.size());
95 labelList oldToNewPoints(origPoints.size(), -1);
104 if (
field[own[facei]] > lowerThreshold)
106 if (
field[nei[facei]] < lowerThreshold)
111 else if (
field[nei[facei]] > lowerThreshold)
117 if (
field[own[facei]] < upperThreshold)
119 if (
field[nei[facei]] > upperThreshold)
124 else if (
field[nei[facei]] < upperThreshold)
132 const face&
f = origFaces[facei];
134 for (
const label pointi :
f)
136 if (oldToNewPoints[pointi] == -1)
138 oldToNewPoints[pointi] =
nPoints++;
153 surfFace =
f.reverseFace();
180 const polyPatch&
p =
bMesh[patchi];
183 zone.start() = nFaces;
187 isA<emptyPolyPatch>(
p)
194 label facei =
p.start();
201 field[own[facei]] > lowerThreshold
202 &&
field[own[facei]] < upperThreshold
205 const face&
f = origFaces[facei];
206 for (
const label pointi :
f)
208 if (oldToNewPoints[pointi] == -1)
210 oldToNewPoints[pointi] =
nPoints++;
214 label
cellId = own[facei];
234 zone.size() =
surfFaces.size() - zone.start();
250 forAll(oldToNewPoints, pointi)
252 if (oldToNewPoints[pointi] >= 0)
254 surfPoints[oldToNewPoints[pointi]] = origPoints[pointi];
274 const scalar lowerThreshold,
275 const scalar upperThreshold,
276 const bool triangulate
281 if (lowerThreshold > upperThreshold)
284 << lowerThreshold <<
" > " << upperThreshold <<
endl;
287 calculate(
field, lowerThreshold, upperThreshold, triangulate);
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
virtual const pointField & points() const
Return raw points.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const surfZoneList & surfZones() const
Const access to the surface zones.
List< face > & storedFaces()
Non-const access to the faces.
virtual label triangulate()
Triangulate in-place, returning the number of triangles added.
thresholdCellFaces(const polyMesh &mesh, const scalarField &field, const scalar lowerThreshold, const scalar upperThreshold, const bool triangulate=false)
Construct from mesh, field and threshold values.
label nInternalFaces() const
Number of internal faces.
pointField & storedPoints()
Non-const access to global points.
label nFaces() const
Number of mesh faces.
static bool & parRun()
Is this a parallel run?
const List< face > & surfFaces() const
Return const access to the faces.
void append(const T &val)
Append an element at the end of the list.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
surfZoneList & storedZones()
Non-const access to the zones.
word name(const complex &c)
Return string representation of complex.
PrimitivePatch< List< face >, const pointField > bMesh
Holder of faceList and points. (v.s. e.g. primitivePatch which references points)
virtual const labelList & faceOwner() const
Return face owner.
void transfer(List< T > &list)
virtual const faceList & faces() const
Return raw faces.
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
List< face > faceList
A List of faces.
void clear()
Clear the list, i.e. set size to zero.
List< surfZone > surfZoneList
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
virtual const labelList & faceNeighbour() const
Return face neighbour.