Go to the documentation of this file.
50 topoSetFaceZoneSource,
57 topoSetFaceZoneSource,
67 planeToFaceZone::typeName,
68 "\n Usage: planeToFaceZone (px py pz) (nx ny nz) include\n\n"
69 " Select faces for which the adjacent cell centres lie on opposite "
78 Foam::planeToFaceZone::faceActionNames_
80 { faceAction::ALL,
"all" },
81 { faceAction::CLOSEST,
"closest" },
87 void Foam::planeToFaceZone::combine(faceZoneSet& fzSet,
const bool add)
const
95 cellIsAbovePlane.set(celli);
109 faceIsOnPlane.set(facei);
117 const label facei =
patch.start() + patchFacei;
120 faceIsOnPlane.set(facei);
127 labelList newSetFaces(faceIsOnPlane.sortedToc());
128 faceIsOnPlane.
clear();
132 if (option_ == faceAction::CLOSEST)
136 labelList newSetFaceRegions(newSetFaces.size(), -1);
142 UIndirectList<face>(
mesh_.
faces(), newSetFaces),
164 procRegionOffset[proci] =
165 procRegionOffset[proci - 1]
166 + procNRegions[proci - 1];
170 for (label& regioni : newSetFaceRegions)
176 nRegions = procRegionOffset.last() + procNRegions.last();
186 forAll(newSetFaces, newSetFacei)
188 const label facei = newSetFaces[newSetFacei];
189 const label regioni = newSetFaceRegions[newSetFacei];
192 meshEdgeRegions[edgei] =
labelList(one{}, regioni);
199 ListOps::appendEqOp<label>(),
205 List<bitSet> regionRegions(nRegions);
206 forAll(newSetFaces, newSetFacei)
208 const label facei = newSetFaces[newSetFacei];
209 const label regioni = newSetFaceRegions[newSetFacei];
213 regionRegions[regioni].set(meshEdgeRegions[edgei]);
215 forAll(regionRegions, regioni)
218 regionRegions[regioni].unset(regioni);
226 forAll(regionRegions, regioni)
228 for (
const label regi : regionRegions[regioni])
231 regionMap[regi] =
min(regionMap[regi], regionMap[regioni]);
241 forAll(regionMap, regioni)
243 if (regionMap[regioni] > regioni0)
246 regionMap[regioni] = regioni0;
251 nRegions = regioni0 + 1;
255 IndirectList<label>(regionMap, newSetFaceRegions);
259 for (
const label regioni : newSetFaceRegions)
261 ++ regionNFaces[regioni];
265 Info<<
" Found " << nRegions <<
" contiguous regions with "
266 << regionNFaces <<
" faces" <<
endl;
270 label selectedRegioni = -1;
275 forAll(newSetFaces, newSetFacei)
277 const label facei = newSetFaces[newSetFacei];
278 const label regioni = newSetFaceRegions[newSetFacei];
283 regionWeights[regioni] += w;
284 regionCentres[regioni] += w*
c;
290 regionCentres /= regionWeights;
301 Info<<
" Selecting region " << selectedRegioni <<
" with "
302 << regionNFaces[selectedRegioni]
303 <<
" faces as the closest to point " << point_ <<
endl;
309 label newSetFacei0 = 0;
310 forAll(newSetFaces, newSetFacei)
312 newSetFaces[newSetFacei0] = newSetFaces[newSetFacei];
314 if (newSetFaceRegions[newSetFacei] == selectedRegioni)
319 newSetFaces.resize(newSetFacei0);
324 DynamicList<label> newAddressing;
325 DynamicList<bool> newFlipMap;
329 newAddressing = fzSet.addressing();
330 newFlipMap = fzSet.flipMap();
333 const auto& exclude = fzSet;
334 for (
const label facei : newSetFaces)
336 if (!exclude.found(facei))
338 newAddressing.append(facei);
346 newAddressing.reserve(fzSet.addressing().size());
347 newFlipMap.reserve(newAddressing.capacity());
350 bitSet exclude(newSetFaces);
351 for (
const label facei : fzSet.addressing())
353 if (!exclude.found(facei))
355 newAddressing.append(facei);
360 fzSet.addressing().transfer(newAddressing);
361 fzSet.flipMap().transfer(newFlipMap);
371 const point& basePoint,
394 faceActionNames_.getOrDefault(
"option",
dict, faceAction::ALL)
407 normal_(checkIs(is)),
408 option_(faceActionNames_.read(checkIs(is)))
420 if (!isA<faceZoneSet>(
set))
423 <<
"Operation only allowed on a faceZoneSet." <<
endl;
433 Info<<
" Adding faces that form a plane at "
434 << point_ <<
" with normal " << normal_ <<
endl;
443 Info<<
" Removing faces that form a plane at "
444 << point_ <<
" with normal " << normal_ <<
endl;
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
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.
faceAction
Enumeration defining the valid options.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
void set(List< bool > &bools, const labelRange &range)
Set the specified range 'on' in a boolList.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
static constexpr const zero Zero
Global zero (0)
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
label nFaces() const
Number of mesh faces.
Class with constructor to add usage string to table.
The topoSetFaceZoneSource is a intermediate class for handling topoSet sources for selecting face zon...
const cellList & cells() const
label nEdges() const
Number of mesh edges.
A topoSetSource to select faces based on the adjacent cell centres spanning a given plane....
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
setAction
Enumeration defining the valid actions.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Create a new set and ADD elements to it.
planeToFaceZone()=delete
No default construct.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
label nCells() const
Number of mesh cells.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
messageStream Info
Information stream (uses stdout - output is on the master only)
AccessType combine(const UList< T > &lists, AccessOp aop=accessOp< T >())
Combines sub-lists into a single list.
virtual const labelList & faceOwner() const
Return face owner.
const labelListList & faceEdges() const
addNamedToRunTimeSelectionTable(topoSetCellSource, badQualityToCell, word, badQuality)
Like faceSet but -reads data from faceZone -updates faceZone when writing.
General set of labels of mesh quantity (points, cells, faces).
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
void add(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
Macros for easy insertion into run-time selection tables.
Subtract elements from the set.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
List< labelList > labelListList
A List of labelList.
const vectorField & cellCentres() const
virtual const faceList & faces() const
Return raw faces.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
virtual void applyToSet(const topoSetSource::setAction action, topoSet &) const
Apply specified action to the topoSet.
const std::string patch
OpenFOAM patch number as a std::string.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const vectorField & faceCentres() const
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
void clear()
Clear the list, i.e. set size to zero.
const dimensionedScalar c
Speed of light in a vacuum.
PrimitivePatch< UIndirectList< face >, const pointField & > uindirectPrimitivePatch
A PrimitivePatch with UIndirectList for the faces, const reference for the point field.
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
vector point
Point is a vector.
const polyMesh & mesh_
Reference to the mesh.
label findMin(const ListType &input, label start=0)
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
virtual const labelList & faceNeighbour() const
Return face neighbour.
static label nProcs(const label communicator=worldComm)
Number of processes in parallel run, and 1 for serial run.
const vectorField & faceAreas() const