polyPatch.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2014 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_(NULL),
96  mePtr_(NULL)
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  readLabel(dict.lookup("nFaces")),
126  readLabel(dict.lookup("startFace"))
127  ),
128  bm.mesh().points()
129  ),
130  start_(readLabel(dict.lookup("startFace"))),
131  boundaryMesh_(bm),
132  faceCellsPtr_(NULL),
133  mePtr_(NULL)
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_(NULL),
167  mePtr_(NULL)
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_(NULL),
194  mePtr_(NULL)
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_(NULL),
221  mePtr_(NULL)
222 {}
223 
224 
226 :
227  patchIdentifier(p),
228  primitivePatch(p),
229  start_(p.start_),
230  boundaryMesh_(p.boundaryMesh_),
231  faceCellsPtr_(NULL),
232  mePtr_(NULL)
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 
296 // Return the patch face neighbour cell centres
298 {
299  tmp<vectorField> tcc(new vectorField(size()));
300  vectorField& cc = tcc();
301 
302  // get reference to global cell centres
303  const vectorField& gcc = boundaryMesh_.mesh().cellCentres();
304 
305  const labelUList& faceCells = this->faceCells();
306 
307  forAll(faceCells, facei)
308  {
309  cc[facei] = gcc[faceCells[facei]];
310  }
311 
312  return tcc;
313 }
314 
315 
317 {
318  if (!faceCellsPtr_)
319  {
320  faceCellsPtr_ = new labelList::subList
321  (
322  patchSlice(boundaryMesh().mesh().faceOwner())
323  );
324  }
325 
326  return *faceCellsPtr_;
327 }
328 
329 
331 {
332  if (!mePtr_)
333  {
334  mePtr_ =
335  new labelList
336  (
338  (
339  boundaryMesh().mesh().edges(),
341  )
342  );
343  }
344 
345  return *mePtr_;
346 }
347 
348 
350 {
353  deleteDemandDrivenData(faceCellsPtr_);
354  deleteDemandDrivenData(mePtr_);
355 }
356 
357 
359 {
360  os.writeKeyword("type") << type() << token::END_STATEMENT << nl;
362  os.writeKeyword("nFaces") << size() << token::END_STATEMENT << nl;
363  os.writeKeyword("startFace") << start() << token::END_STATEMENT << nl;
364 }
365 
366 
368 {}
369 
370 
372 (
374  const primitivePatch&,
375  labelList& faceMap,
376  labelList& rotation
377 ) const
378 {
379  // Nothing changed.
380  return false;
381 }
382 
383 
384 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
385 
387 {
388  clearAddressing();
389 
390  patchIdentifier::operator=(p);
392  start_ = p.start_;
393 }
394 
395 
396 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
397 
399 {
400  p.write(os);
401  os.check("Ostream& operator<<(Ostream& os, const polyPatch& p");
402  return os;
403 }
404 
405 
406 // ************************************************************************* //
int debugSwitch(const char *name, const int defaultValue=0)
Lookup debug switch or add default value.
Definition: debug.C:165
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
Definition: polyPatch.C:367
Pre-declare related SubField type.
Definition: Field.H:61
Buffers for inter-processor communications streams (UOPstream, UIPstream).
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
labelList meshEdges(const edgeList &allEdges, const labelListList &cellEdges, const labelList &faceCells) const
Return labels of patch edges in the global edge list using.
const vectorField::subField faceAreas() const
Return face normals.
Definition: polyPatch.C:290
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
Definition: polyPatch.C:358
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:316
void deleteDemandDrivenData(DataPtr &dataPtr)
A class for handling words, derived from string.
Definition: word.H:59
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
Identifies patch by name, patch index and physical type.
const vectorField & cellCentres() const
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Definition: polyPatch.C:62
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:278
dynamicFvMesh & mesh
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Ostream & operator<<(Ostream &, const edgeMesh &)
Definition: edgeMeshIO.C:133
Namespace for OpenFOAM.
void write(Ostream &) const
Write patchIdentifier as a dictionary.
label readLabel(Istream &is)
Definition: label.H:64
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
const vectorField::subField faceCentres() const
Return face centres.
Definition: polyPatch.C:284
Abstract base class for point-mesh patch fields.
static const char nl
Definition: Ostream.H:260
void setSize(const label)
Reset size of List.
Definition: List.C:318
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
static wordList constraintTypes()
Return a list of all the constraint patch types.
Definition: polyPatch.C:252
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1035
Foam::polyBoundaryMesh.
volScalarField & p
Definition: createFields.H:51
static int disallowGenericPolyPatch
Debug switch to disallow the use of genericPolyPatch.
Definition: polyPatch.H:131
tmp< vectorField > faceCellCentres() const
Return face cell centres.
Definition: polyPatch.C:297
Field< vector > vectorField
Specialisation of Field<T> for vector.
#define forAll(list, i)
Definition: UList.H:421
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
static bool constraintType(const word &pt)
Return true if the given type is a constraint type.
Definition: polyPatch.C:246
virtual ~polyPatch()
Destructor.
Definition: polyPatch.C:238
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
Macros for easy insertion into run-time selection tables.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
A list of faces which address into the list of points.
const List< T >::subList patchSlice(const UList< T > &l) const
Slice list to patch.
Definition: polyPatch.H:333
virtual void clearGeom()
Clear geometry.
Definition: polyPatch.C:69
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
void operator=(const polyPatch &)
Assignment.
Definition: polyPatch.C:386
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
const labelListList & pointEdges() const
Return point-edge addressing.
const polyMesh & mesh() const
Return the mesh reference.
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
SubList< label > subList
Declare type of subList.
Definition: List.H:153
A List obtained as a section of another List.
Definition: SubList.H:53
PrimitivePatch< face, SubList, const pointField & > primitivePatch
Addressing for a faceList slice.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1060
void operator=(const PrimitivePatch< Face, FaceList, PointField, PointType > &)
Assignment.
SubList< face > faceSubList
Definition: faceListFwd.H:44
virtual void movePoints(const Field< PointType > &)
Correct patch after moving points.
virtual void clearAddressing()
Clear addressing.
Definition: polyPatch.C:349
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
A class for managing temporary objects.
Definition: PtrList.H:118
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
Definition: polyPatch.C:372
static const word null
An empty word.
Definition: word.H:77
defineTypeNameAndDebug(combustionModel, 0)
virtual void movePoints(PstreamBuffers &, const pointField &p)
Correct patches after moving points.
Definition: polyPatch.C:57
const labelList & meshEdges() const
Return global edge index for local edges.
Definition: polyPatch.C:330