mappedPolyPatch.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2012 OpenFOAM Foundation
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 Class
27  Foam::mappedPolyPatch
28 
29 Description
30  Determines a mapping between patch face centres and mesh cell or face
31  centres and processors they're on.
32 
33 Note
34  Storage is not optimal. It stores all face centres and cells on all
35  processors to keep the addressing calculation simple.
36 
37 SourceFiles
38  mappedPolyPatch.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef mappedPolyPatch_H
43 #define mappedPolyPatch_H
44 
45 #include "polyPatch.H"
46 #include "mappedPatchBase.H"
47 
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 
54 class polyMesh;
55 
56 /*---------------------------------------------------------------------------*\
57  Class mappedPolyPatch Declaration
58 \*---------------------------------------------------------------------------*/
59 
60 class mappedPolyPatch
61 :
62  public polyPatch,
63  public mappedPatchBase
64 {
65 
66 protected:
67 
68  //- Initialise the calculation of the patch geometry
69  virtual void initGeometry(PstreamBuffers&);
70 
71  //- Calculate the patch geometry
72  virtual void calcGeometry(PstreamBuffers&);
73 
74  //- Initialise the patches for moving points
75  virtual void initMovePoints(PstreamBuffers&, const pointField&);
76 
77  //- Correct patches after moving points
78  virtual void movePoints(PstreamBuffers&, const pointField&);
79 
80  //- Initialise the update of the patch topology
81  virtual void initUpdateMesh(PstreamBuffers&);
82 
83  //- Update of the patch topology
84  virtual void updateMesh(PstreamBuffers&);
85 
86 
87 public:
88 
89  //- Runtime type information
90  TypeName("mappedPatch");
91 
92 
93  // Constructors
94 
95  //- Construct from components
97  (
98  const word& name,
99  const label size,
100  const label start,
101  const label index,
102  const polyBoundaryMesh& bm,
103  const word& patchType
104  );
105 
106  //- Construct from components
108  (
109  const word& name,
110  const label size,
111  const label start,
112  const label index,
113  const word& sampleRegion,
115  const word& samplePatch,
116  const vectorField& offset,
117  const polyBoundaryMesh& bm
118  );
119 
120  //- Construct from components. Uniform offset.
122  (
123  const word& name,
124  const label size,
125  const label start,
126  const label index,
127  const word& sampleRegion,
129  const word& samplePatch,
130  const vector& offset,
131  const polyBoundaryMesh& bm
132  );
133 
134  //- Construct from dictionary
136  (
137  const word& name,
138  const dictionary& dict,
139  const label index,
140  const polyBoundaryMesh& bm,
141  const word& patchType
142  );
143 
144  //- Construct as copy, resetting the boundary mesh
146  (
147  const mappedPolyPatch&,
148  const polyBoundaryMesh&
149  );
150 
151  //- Construct given the original patch and resetting the
152  // face list and boundary mesh information
154  (
155  const mappedPolyPatch& pp,
156  const polyBoundaryMesh& bm,
157  const label index,
158  const label newSize,
159  const label newStart
160  );
161 
162  //- Construct given the original patch and a map
164  (
165  const mappedPolyPatch& pp,
166  const polyBoundaryMesh& bm,
167  const label index,
168  const labelUList& mapAddressing,
169  const label newStart
170  );
171 
172  //- Construct and return a clone, resetting the boundary mesh
173  virtual autoPtr<polyPatch> clone(const polyBoundaryMesh& bm) const
174  {
175  return autoPtr<polyPatch>(new mappedPolyPatch(*this, bm));
176  }
177 
178  //- Construct and return a clone, resetting the face list
179  // and boundary mesh
180  virtual autoPtr<polyPatch> clone
181  (
182  const polyBoundaryMesh& bm,
183  const label index,
184  const label newSize,
185  const label newStart
186  ) const
187  {
188  return autoPtr<polyPatch>
189  (
190  new mappedPolyPatch(*this, bm, index, newSize, newStart)
191  );
192  }
193 
194  //- Construct and return a clone, resetting the face list
195  // and boundary mesh
196  virtual autoPtr<polyPatch> clone
197  (
198  const polyBoundaryMesh& bm,
199  const label index,
200  const labelUList& mapAddressing,
201  const label newStart
202  ) const
203  {
204  return autoPtr<polyPatch>
205  (
206  new mappedPolyPatch
207  (
208  *this,
209  bm,
211  mapAddressing,
212  newStart
213  )
214  );
215  }
216 
217  // Implicit treatment functions
218 
219  //- Return number of new internal of this polyPatch faces
220  virtual void newInternalProcFaces
221  (
222  label& iFaces,
223  label& pFaces
224  ) const
225  {
226  label nbrFaces =
228  iFaces = patch_.size();
229  pFaces = max(iFaces - nbrFaces, 0);
230  }
231 
232 
233  //- Return nbrCells
234  virtual const labelUList& nbrCells() const
235  {
237  }
238 
239  //- Return nbr patch ID
240  virtual label neighbPolyPatchID() const
241  {
243  }
244 
245  //- Return collocated faces map
246  virtual refPtr<labelListList> mapCollocatedFaces() const
247  {
248  refPtr<labelListList> tMap(new labelListList(patch_.size()));
249  labelListList& map = tMap.ref();
250  forAll(map, i)
251  {
252  labelList& subMap = map[i];
253  subMap.setSize(1);
254  subMap[0] = i;
255  }
256  return tMap;
257  }
258 
259  //- Return implicit master
260  virtual bool masterImplicit() const
261  {
262  return owner();
263  }
264 
265  //- Return neigh region ID
266  virtual word neighbRegionID() const
267  {
268  return sampleRegion_;
269  }
271 
272  //- Destructor
273  virtual ~mappedPolyPatch();
274 
275 
276  // Member functions
277 
278  //- Write the polyPatch data as a dictionary
279  virtual void write(Ostream&) const;
280 };
281 
282 
283 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
284 
285 } // End namespace Foam
287 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
288 
289 #endif
290 
291 // ************************************************************************* //
word samplePatch_
Patch (if in sampleMode NEARESTPATCH*)
virtual label neighbPolyPatchID() const
Return nbr patch ID.
dictionary dict
virtual const labelUList & nbrCells() const
Return nbrCells.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
const vector & offset() const noexcept
Offset vector (from patch faces to destination mesh objects)
word sampleRegion_
Region to sample.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:129
label start() const noexcept
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:441
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
const polyPatch & patch_
Patch to sample.
const polyPatch & lookupPatch(const word &sampleRegion, const word &samplePatch) const
Lookup patch.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:40
mappedPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType)
Construct from components.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
const mapDistribute & map() const
Return reference to the parallel distribution map.
List< labelList > labelListList
List of labelList.
Definition: labelList.H:38
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:421
virtual word neighbRegionID() const
Return neigh region ID.
TypeName("mappedPatch")
Runtime type information.
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:365
A class for handling words, derived from Foam::string.
Definition: word.H:63
virtual void movePoints(PstreamBuffers &, const pointField &)
Correct patches after moving points.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO...
const word & name() const noexcept
The patch name.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:56
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
virtual bool masterImplicit() const
Return implicit master.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=cellModel::ref(cellModel::HEX);labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){ pointMap[i]=i;}for(label i=0;i< nPoints;i++){ if(f[i] > 0.0) { hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells].reset(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei) { if(edges[ei].mag(points)< SMALL) { label start=pointMap[edges[ei].start()];while(start !=pointMap[start]) { start=pointMap[start];} label end=pointMap[edges[ei].end()];while(end !=pointMap[end]) { end=pointMap[end];} label minLabel=min(start, end);pointMap[start]=pointMap[end]=minLabel;} } cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){ cellShape &cs=cellShapes[celli];forAll(cs, i) { cs[i]=pointMap[cs[i]];} cs.collapse();}label bcIDs[11]={-1, 0, 2, 4, -1, 5, -1, 6, 7, 8, 9};const label nBCs=12;const word *kivaPatchTypes[nBCs]={ &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &symmetryPolyPatch::typeName, &wedgePolyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &symmetryPolyPatch::typeName, &oldCyclicPolyPatch::typeName};enum patchTypeNames{ PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={ "piston", "valve", "liner", "cylinderHead", "axis", "wedge", "inflow", "outflow", "presin", "presout", "symmetryPlane", "cyclic"};List< SLList< face > > pFaces[nBCs]
Definition: readKivaGrid.H:235
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm) const
Construct and return a clone, resetting the boundary mesh.
sampleMode mode() const noexcept
What to sample.
virtual refPtr< labelListList > mapCollocatedFaces() const
Return collocated faces map.
bool owner() const
Is it owner.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
sampleMode
Mesh items to sample.
const word & samplePatch() const
Patch (only if NEARESTPATCHFACE)
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
Definition: HashPtrTable.H:48
label index() const noexcept
The index of this patch in the boundaryMesh.
List< label > labelList
A List of labels.
Definition: List.H:62
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:69
const word & sampleRegion() const
Region to sample.
virtual ~mappedPolyPatch()
Destructor.
virtual void newInternalProcFaces(label &iFaces, label &pFaces) const
Return number of new internal of this polyPatch faces.
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
Namespace for OpenFOAM.