polyPatch.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "polyPatch.H"
28 #include "polyBoundaryMesh.H"
29 #include "polyMesh.H"
30 #include "primitiveMesh.H"
31 #include "SubField.H"
32 #include "entry.H"
33 #include "dictionary.H"
34 #include "pointPatchField.H"
35 
36 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37 
38 namespace Foam
39 {
40  defineTypeNameAndDebug(polyPatch, 0);
41 
43  (
44  debug::debugSwitch("disallowGenericPolyPatch", 0)
45  );
46 
47  defineRunTimeSelectionTable(polyPatch, word);
48  defineRunTimeSelectionTable(polyPatch, dictionary);
49 
50  addToRunTimeSelectionTable(polyPatch, polyPatch, word);
51  addToRunTimeSelectionTable(polyPatch, polyPatch, dictionary);
52 }
53 
54 
55 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
56 
58 {
60 }
61 
63 {
65  clearAddressing();
66 }
67 
68 
70 {
72 }
73 
74 
75 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
76 
78 (
79  const word& name,
80  const label size,
81  const label start,
82  const label index,
83  const polyBoundaryMesh& bm,
84  const word& patchType
85 )
86 :
87  patchIdentifier(name, index),
89  (
90  faceSubList(bm.mesh().faces(), size, start),
91  bm.mesh().points()
92  ),
93  start_(start),
94  boundaryMesh_(bm),
95  faceCellsPtr_(nullptr),
96  mePtr_(nullptr)
97 {
98  if
99  (
100  patchType != word::null
101  && constraintType(patchType)
102  && findIndex(inGroups(), patchType) == -1
103  )
104  {
105  inGroups().append(patchType);
106  }
107 }
108 
109 
111 (
112  const word& name,
113  const dictionary& dict,
114  const label index,
115  const polyBoundaryMesh& bm,
116  const word& patchType
117 )
118 :
119  patchIdentifier(name, dict, index),
121  (
123  (
124  bm.mesh().faces(),
125  dict.lookup<label>("nFaces"),
126  dict.lookup<label>("startFace")
127  ),
128  bm.mesh().points()
129  ),
130  start_(dict.lookup<label>("startFace")),
131  boundaryMesh_(bm),
132  faceCellsPtr_(nullptr),
133  mePtr_(nullptr)
134 {
135  if
136  (
137  patchType != word::null
138  && constraintType(patchType)
139  && findIndex(inGroups(), patchType) == -1
140  )
141  {
142  inGroups().append(patchType);
143  }
144 }
145 
146 
148 (
149  const polyPatch& pp,
150  const polyBoundaryMesh& bm
151 )
152 :
153  patchIdentifier(pp),
155  (
157  (
158  bm.mesh().faces(),
159  pp.size(),
160  pp.start()
161  ),
162  bm.mesh().points()
163  ),
164  start_(pp.start()),
165  boundaryMesh_(bm),
166  faceCellsPtr_(nullptr),
167  mePtr_(nullptr)
168 {}
169 
170 
172 (
173  const polyPatch& pp,
174  const polyBoundaryMesh& bm,
175  const label index,
176  const label newSize,
177  const label newStart
178 )
179 :
180  patchIdentifier(pp, index),
182  (
184  (
185  bm.mesh().faces(),
186  newSize,
187  newStart
188  ),
189  bm.mesh().points()
190  ),
191  start_(newStart),
192  boundaryMesh_(bm),
193  faceCellsPtr_(nullptr),
194  mePtr_(nullptr)
195 {}
196 
197 
199 (
200  const polyPatch& pp,
201  const polyBoundaryMesh& bm,
202  const label index,
203  const labelUList& mapAddressing,
204  const label newStart
205 )
206 :
207  patchIdentifier(pp, index),
209  (
211  (
212  bm.mesh().faces(),
213  mapAddressing.size(),
214  newStart
215  ),
216  bm.mesh().points()
217  ),
218  start_(newStart),
219  boundaryMesh_(bm),
220  faceCellsPtr_(nullptr),
221  mePtr_(nullptr)
222 {}
223 
224 
226 :
227  patchIdentifier(p),
228  primitivePatch(p),
229  start_(p.start_),
230  boundaryMesh_(p.boundaryMesh_),
231  faceCellsPtr_(nullptr),
232  mePtr_(nullptr)
233 {}
234 
235 
236 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
237 
239 {
240  clearAddressing();
241 }
242 
243 
244 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
245 
247 {
249 }
250 
251 
253 {
254  wordList cTypes(dictionaryConstructorTablePtr_->size());
255 
256  label i = 0;
257 
258  for
259  (
260  dictionaryConstructorTable::iterator cstrIter =
261  dictionaryConstructorTablePtr_->begin();
262  cstrIter != dictionaryConstructorTablePtr_->end();
263  ++cstrIter
264  )
265  {
266  if (constraintType(cstrIter.key()))
267  {
268  cTypes[i++] = cstrIter.key();
269  }
270  }
271 
272  cTypes.setSize(i);
273 
274  return cTypes;
275 }
276 
277 
279 {
280  return boundaryMesh_;
281 }
282 
283 
285 {
286  return patchSlice(boundaryMesh().mesh().faceCentres());
287 }
288 
289 
291 {
292  return patchSlice(boundaryMesh().mesh().faceAreas());
293 }
294 
295 
297 {
299 }
300 
301 
303 {
304  tmp<vectorField> tcc(new vectorField(size()));
305  vectorField& cc = tcc.ref();
306 
307  // get reference to global cell centres
308  const vectorField& gcc = boundaryMesh_.mesh().cellCentres();
309 
310  const labelUList& faceCells = this->faceCells();
311 
312  forAll(faceCells, facei)
313  {
314  cc[facei] = gcc[faceCells[facei]];
315  }
316 
317  return tcc;
318 }
319 
320 
322 {
323  if (!faceCellsPtr_)
324  {
325  faceCellsPtr_ = new labelList::subList
326  (
327  patchSlice(boundaryMesh().mesh().faceOwner())
328  );
329  }
330 
331  return *faceCellsPtr_;
332 }
333 
334 
336 {
337  if (!mePtr_)
338  {
339  mePtr_ =
340  new labelList
341  (
343  (
344  boundaryMesh().mesh().edges(),
346  )
347  );
348  }
349 
350  return *mePtr_;
351 }
352 
353 
355 {
358  deleteDemandDrivenData(faceCellsPtr_);
359  deleteDemandDrivenData(mePtr_);
360 }
361 
362 
364 {
365  writeEntry(os, "type", type());
367  writeEntry(os, "nFaces", size());
368  writeEntry(os, "startFace", start());
369 }
370 
371 
373 {}
374 
375 
377 (
379  const primitivePatch&,
380  labelList& faceMap,
381  labelList& rotation
382 ) const
383 {
384  // Nothing changed.
385  return false;
386 }
387 
388 
389 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
390 
392 {
393  clearAddressing();
394 
395  patchIdentifier::operator=(p);
397  start_ = p.start_;
398 }
399 
400 
401 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
402 
404 {
405  p.write(os);
406  os.check("Ostream& operator<<(Ostream& os, const polyPatch& p");
407  return os;
408 }
409 
410 
411 // ************************************************************************* //
const labelListList & pointEdges() const
Return point-edge addressing.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual void movePoints(const Field< PointType > &)
Correct patch after moving points.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
static int disallowGenericPolyPatch
Debug switch to disallow the use of genericPolyPatch.
Definition: polyPatch.H:131
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
Identifies patch by name, patch index and physical type.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:278
SubList< label > subList
Declare type of subList.
Definition: List.H:199
virtual void clearAddressing()
Clear addressing.
Definition: polyPatch.C:354
Pre-declare related SubField type.
Definition: Field.H:60
Macros for easy insertion into run-time selection tables.
virtual void clearGeom()
Clear geometry.
Definition: polyPatch.C:69
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
Definition: polyPatch.C:363
virtual void movePoints(PstreamBuffers &, const pointField &p)
Correct patches after moving points.
Definition: polyPatch.C:57
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
static bool constraintType(const word &pt)
Return true if the given type is a constraint type.
Definition: polyPatch.C:246
Abstract base class for point-mesh patch fields.
A list of faces which address into the list of points.
A List obtained as a section of another List.
Definition: SubList.H:53
dynamicFvMesh & mesh
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:321
A class for handling words, derived from string.
Definition: word.H:59
virtual ~polyPatch()
Destructor.
Definition: polyPatch.C:238
int debugSwitch(const char *name, const int defaultValue=0)
Lookup debug switch or add default value.
Definition: debug.C:178
static const word null
An empty word.
Definition: word.H:77
const polyMesh & mesh() const
Return the mesh reference.
tmp< vectorField > faceCellCentres() const
Return face cell centres.
Definition: polyPatch.C:302
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
List< label > labelList
A List of labels.
Definition: labelList.H:56
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Definition: polyPatch.C:62
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
polyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType)
Construct from components.
Definition: polyPatch.C:78
const vectorField & cellCentres() const
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
void operator=(const polyPatch &)
Assignment.
Definition: polyPatch.C:391
Foam::polyBoundaryMesh.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
const vectorField::subField faceAreas() const
Return face areas.
Definition: polyPatch.C:290
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
defineTypeNameAndDebug(combustionModel, 0)
const labelList & meshEdges() const
Return global edge index for local edges.
Definition: polyPatch.C:335
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
Addressing for a faceList slice.
void write(Ostream &) const
Write patchIdentifier as a dictionary.
static wordList constraintTypes()
Return a list of all the constraint patch types.
Definition: polyPatch.C:252
void setSize(const label)
Reset size of List.
Definition: List.C:281
const scalarField::subField magFaceAreas() const
Return face area magnitudes.
Definition: polyPatch.C:296
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
Definition: polyPatch.C:377
labelList meshEdges(const edgeList &allEdges, const labelListList &cellEdges, const labelList &faceCells) const
Return labels of patch edges in the global edge list using.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:303
Ostream & operator<<(Ostream &, const ensightPart &)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
SubList< face > faceSubList
Definition: faceListFwd.H:44
Field< vector > vectorField
Specialisation of Field<T> for vector.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
Definition: polyPatch.C:372
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
void deleteDemandDrivenData(DataPtr &dataPtr)
void operator=(const PrimitivePatch< FaceList, PointField > &)
Assignment operator.
const List< T >::subList patchSlice(const UList< T > &l) const
Slice list to patch.
Definition: polyPatch.H:336
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812
const vectorField::subField faceCentres() const
Return face centres.
Definition: polyPatch.C:284