Loading...
Searching...
No Matches
mappedWallPolyPatch.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-------------------------------------------------------------------------------
10License
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
26Class
27 Foam::mappedWallPolyPatch
28
29Description
30 Determines a mapping between patch face centres and mesh cell or face
31 centres and processors they're on.
32
33Note
34 Storage is not optimal. It stores all face centres and cells on all
35 processors to keep the addressing calculation simple.
36
37SourceFiles
38 mappedWallPolyPatch.C
39
40\*---------------------------------------------------------------------------*/
41
42#ifndef mappedWallPolyPatch_H
43#define mappedWallPolyPatch_H
44
45#include "wallPolyPatch.H"
46#include "mappedPatchBase.H"
47
48
49// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50
51namespace Foam
52{
53
54class polyMesh;
56/*---------------------------------------------------------------------------*\
57 Class mappedWallPolyPatch Declaration
58\*---------------------------------------------------------------------------*/
59
61:
62 public wallPolyPatch,
63 public mappedPatchBase
64{
65
66protected:
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
87public:
88
89 //- Runtime type information
90 TypeName("mappedWall");
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 mappedWallPolyPatch&,
148 const polyBoundaryMesh&
149 );
150
151 //- Construct given the original patch and resetting the
152 //- face list and boundary mesh information
154 (
155 const mappedWallPolyPatch& 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 mappedWallPolyPatch& 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 mappedWallPolyPatch(*this, bm));
176 }
177
178 //- Construct and return a clone, resetting the face list
179 //- and boundary mesh
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 (
191 (
192 *this,
193 bm,
194 index,
195 newSize,
196 newStart
197 )
198 );
199 }
200
201 //- Construct and return a clone, resetting the face list
202 //- and boundary mesh
204 (
205 const polyBoundaryMesh& bm,
206 const label index,
207 const labelUList& mapAddressing,
208 const label newStart
209 ) const
210 {
211 return autoPtr<polyPatch>
212 (
214 (
215 *this,
216 bm,
217 index,
218 mapAddressing,
219 newStart
220 )
221 );
222 }
223
224 // Implicit treatment functions
225
226 //- Return number of new internal of this polyPatch faces
227 virtual void newInternalProcFaces(label& iFaces, label& pFaces) const
228 {
229 label nbrFaces = lookupPatch(sampleRegion_, samplePatch_).size();
230 iFaces = patch_.size();
231 pFaces = max(iFaces - nbrFaces, 0);
233
234
235 //- Return nbrCells
236 virtual const labelUList& nbrCells() const
237 {
239 }
240
241 //- Return nbr patch ID
242 virtual label neighbPolyPatchID() const
243 {
245 }
246
247 //- Return collocated faces map
248 virtual refPtr<labelListList> mapCollocatedFaces() const
249 {
250 refPtr<labelListList> tMap(new labelListList(patch_.size()));
251 labelListList& map = tMap.ref();
252 forAll (map, i)
253 {
254 labelList& subMap = map[i];
255 subMap.setSize(1);
256 subMap[0] = i;
257 }
258 return tMap;
259 }
260
261 //- Return implicit master
262 virtual bool masterImplicit() const
263 {
264 return owner();
265 }
266
267 //- Return neigh region ID
268 virtual word neighbRegionID() const
270 return sampleRegion_;
271 }
272
273
274 //- Destructor
275 virtual ~mappedWallPolyPatch();
276
278 // Member functions
279
280 //- Write the polyPatch data as a dictionary
281 virtual void write(Ostream&) const;
282};
283
284
285// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
286
287} // End namespace Foam
288
289// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
290
291#endif
292
293// ************************************************************************* //
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
void setSize(label n)
Alias for resize().
Definition List.H:536
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition Ostream.H:59
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void size(const label n)
Definition UList.H:118
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition autoPtr.H:65
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition dictionary.H:133
const mapDistribute & map() const
Return reference to the parallel distribution map.
const vector & offset() const noexcept
Offset vector (from patch faces to destination mesh objects).
bool owner() const
Is it owner.
const polyPatch & patch_
Patch to sample.
sampleMode
Mesh items to sample.
const polyPatch & lookupPatch(const word &sampleRegion, const word &samplePatch) const
Lookup patch.
word sampleRegion_
Region to sample.
word samplePatch_
Patch (if in sampleMode NEARESTPATCH*).
mappedPatchBase(const polyPatch &)
Construct from patch.
const word & samplePatch() const
Patch (only if NEARESTPATCHFACE).
sampleMode mode() const noexcept
What to sample.
const word & sampleRegion() const
Region to sample.
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm) const
Construct and return a clone, resetting the boundary mesh.
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
mappedWallPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType)
Construct from components.
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
virtual label neighbPolyPatchID() const
Return nbr patch ID.
virtual bool masterImplicit() const
Return implicit master.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm, const label index, const label newSize, const label newStart) const
Construct and return a clone, resetting the face list and boundary mesh.
virtual const labelUList & nbrCells() const
Return nbrCells.
TypeName("mappedWall")
Runtime type information.
virtual refPtr< labelListList > mapCollocatedFaces() const
Return collocated faces map.
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm, const label index, const labelUList &mapAddressing, const label newStart) const
Construct and return a clone, resetting the face list and boundary mesh.
virtual ~mappedWallPolyPatch()
Destructor.
virtual void movePoints(PstreamBuffers &, const pointField &)
Correct patches after moving points.
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
virtual void newInternalProcFaces(label &iFaces, label &pFaces) const
Return number of new internal of this polyPatch faces.
virtual word neighbRegionID() const
Return neigh region ID.
label index() const noexcept
The index of this patch in the boundaryMesh.
A polyBoundaryMesh is a polyPatch list with registered IO, a reference to the associated polyMesh,...
Mesh consisting of general polyhedral cells.
Definition polyMesh.H:79
label start() const noexcept
Return start label of this patch in the polyMesh face list.
Definition polyPatch.H:446
const labelUList & faceCells() const
Return face-cell addressing.
Definition polyPatch.C:401
A class for managing references or pointers (no reference counting).
Definition refPtr.H:54
T & ref() const
Return non-const reference to the contents of a non-null managed pointer.
Definition refPtrI.H:230
wallPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType)
Construct from components.
A class for handling words, derived from Foam::string.
Definition word.H:66
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition hashSets.C:40
List< labelList > labelListList
List of labelList.
Definition labelList.H:38
List< label > labelList
A List of labels.
Definition List.H:62
Field< vector > vectorField
Specialisation of Field<T> for vector.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
Definition exprTraits.C:127
vectorField pointField
pointField is a vectorField.
Vector< scalar > vector
Definition vector.H:57
UList< label > labelUList
A UList of labels.
Definition UList.H:75
runTime write()
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];} pointMap[start]=pointMap[end]=Foam::min(start, end);} } 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};constexpr 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< DynamicList< face > > pFaces[nBCs]
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition stdFoam.H:299
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition typeInfo.H:68