extrudePatchMesh.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 "extrudePatchMesh.H"
27 
28 #include "createShellMesh.H"
29 #include "polyTopoChange.H"
30 #include "wallPolyPatch.H"
31 #include "emptyPolyPatch.H"
32 #include "wedgePolyPatch.H"
33 
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 defineTypeNameAndDebug(extrudePatchMesh, 0);
43 
44 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45 
47 (
48  const fvMesh& mesh,
49  const fvPatch& patch,
50  const dictionary& dict,
51  const word regionName,
52  const List<polyPatch*>& regionPatches
53 )
54 :
55  fvMesh
56  (
57  IOobject
58  (
59  regionName,
60  mesh.facesInstance(),
61  mesh,
64  true
65  ),
66  pointField(),
67  faceList(),
68  labelList(),
69  labelList(),
70  false
71  ),
72  extrudedPatch_(patch.patch()),
73  dict_(dict)
74 {
75  extrudeMesh(regionPatches);
76 }
77 
78 
80 (
81  const fvMesh& mesh,
82  const fvPatch& patch,
83  const dictionary& dict,
84  const word regionName
85 )
86 :
87  fvMesh
88  (
89  IOobject
90  (
91  regionName,
92  mesh.facesInstance(),
93  mesh,
96  true
97  ),
98  pointField(),
99  faceList(),
100  labelList(),
101  labelList(),
102  false
103  ),
104  extrudedPatch_(patch.patch()),
105  dict_(dict)
106 {
107 
108  List<polyPatch*> regionPatches(3);
109  List<word> patchNames(regionPatches.size());
110  List<word> patchTypes(regionPatches.size());
111  PtrList<dictionary> dicts(regionPatches.size());
112 
113  forAll(dicts, patchi)
114  {
115  if (!dicts.set(patchi))
116  {
117  dicts.set(patchi, new dictionary());
118  }
119  }
120 
121  dicts[bottomPatchID] = dict_.subDict("bottomCoeffs");
122  dicts[sidePatchID] = dict_.subDict("sideCoeffs");
123  dicts[topPatchID] = dict_.subDict("topCoeffs");
124 
125  forAll(dicts, patchi)
126  {
127  dicts[patchi].lookup("name") >> patchNames[patchi];
128  dicts[patchi].lookup("type") >> patchTypes[patchi];
129  }
130 
131  forAll(regionPatches, patchi)
132  {
133  dictionary& patchDict = dicts[patchi];
134  patchDict.set("nFaces", 0);
135  patchDict.set("startFace", 0);
136 
137  regionPatches[patchi] = polyPatch::New
138  (
140  patchDict,
141  patchi,
142  mesh.boundaryMesh()
143  ).ptr();
144 
145  }
146 
147  extrudeMesh(regionPatches);
148 
149 }
150 
151 
152 void extrudePatchMesh::extrudeMesh(const List<polyPatch*>& regionPatches)
153 {
154  if (this->boundaryMesh().size() == 0)
155  {
156  bool columnCells = readBool(dict_.lookup("columnCells"));
157 
158  PackedBoolList nonManifoldEdge(extrudedPatch_.nEdges());
159  for (label edgeI = 0; edgeI < extrudedPatch_.nInternalEdges(); edgeI++)
160  {
161  if (columnCells)
162  {
163  nonManifoldEdge[edgeI] = true;
164  }
165  }
166 
168 
169  faceList pointGlobalRegions;
170  faceList pointLocalRegions;
171  labelList localToGlobalRegion;
172 
173  const primitiveFacePatch pp
174  (
175  extrudedPatch_, extrudedPatch_.points()
176  );
177 
179  (
180  this->globalData(),
181  pp,
182  nonManifoldEdge,
183  false,
184 
185  pointGlobalRegions,
186  pointLocalRegions,
187  localToGlobalRegion
188  );
189 
190 
191  // Per local region an originating point
192  labelList localRegionPoints(localToGlobalRegion.size());
193  forAll(pointLocalRegions, facei)
194  {
195  const face& f = extrudedPatch_.localFaces()[facei];
196  const face& pRegions = pointLocalRegions[facei];
197  forAll(pRegions, fp)
198  {
199  localRegionPoints[pRegions[fp]] = f[fp];
200  }
201  }
202 
203  // Calculate region normals by reducing local region normals
204  pointField localRegionNormals(localToGlobalRegion.size());
205  {
206  pointField localSum(localToGlobalRegion.size(), Zero);
207 
208  forAll(pointLocalRegions, facei)
209  {
210  const face& pRegions = pointLocalRegions[facei];
211  forAll(pRegions, fp)
212  {
213  label localRegionI = pRegions[fp];
214  localSum[localRegionI] +=
215  extrudedPatch_.faceNormals()[facei];
216  }
217  }
218 
219  Map<point> globalSum(2*localToGlobalRegion.size());
220 
221  forAll(localSum, localRegionI)
222  {
223  label globalRegionI = localToGlobalRegion[localRegionI];
224  globalSum.insert(globalRegionI, localSum[localRegionI]);
225  }
226 
227  // Reduce
229  Pstream::mapCombineScatter(globalSum);
230 
231  forAll(localToGlobalRegion, localRegionI)
232  {
233  label globalRegionI = localToGlobalRegion[localRegionI];
234  localRegionNormals[localRegionI] = globalSum[globalRegionI];
235  }
236  localRegionNormals /= mag(localRegionNormals);
237  }
238 
239 
240  // Per local region an extrusion direction
241  vectorField firstDisp(localToGlobalRegion.size());
242  forAll(firstDisp, regionI)
243  {
244  // const point& regionPt = regionCentres[regionI];
245  const point& regionPt = extrudedPatch_.points()
246  [
247  extrudedPatch_.meshPoints()
248  [
249  localRegionPoints[regionI]
250  ]
251  ];
252  const vector& n = localRegionNormals[regionI];
253  firstDisp[regionI] = model_()(regionPt, n, 1) - regionPt;
254  }
255 
256 
257  // Extrude engine
258  createShellMesh extruder
259  (
260  pp,
261  pointLocalRegions,
262  localRegionPoints
263  );
264 /*
265  List<polyPatch*> regionPatches(3);
266  List<word> patchNames(regionPatches.size());
267  List<word> patchTypes(regionPatches.size());
268  PtrList<dictionary> dicts(regionPatches.size());
269 
270  forAll(dicts, patchi)
271  {
272  if (!dicts.set(patchi))
273  {
274  dicts.set(patchi, new dictionary());
275  }
276  }
277 
278  dicts[bottomPatchID] = dict_.subDict("bottomCoeffs");
279  dicts[sidePatchID] = dict_.subDict("sideCoeffs");
280  dicts[topPatchID] = dict_.subDict("topCoeffs");
281 
282  forAll(dicts, patchi)
283  {
284  dicts[patchi].lookup("name") >> patchNames[patchi];
285  dicts[patchi].lookup("type") >> patchTypes[patchi];
286  }
287 
288  forAll(regionPatches, patchi)
289  {
290  dictionary& patchDict = dicts[patchi];
291  patchDict.set("nFaces", 0);
292  patchDict.set("startFace", 0);
293 
294  regionPatches[patchi] = polyPatch::New
295  (
296  patchNames[patchi],
297  patchDict,
298  patchi,
299  mesh.boundaryMesh()
300  ).ptr();
301 
302  }
303 */
304  this->clearOut();
305  this->removeFvBoundary();
306  this->addFvPatches(regionPatches, true);
307 
308 
309  // At this point we have a valid mesh with 3 patches and zero cells.
310  // Determine:
311  // - per face the top and bottom patch (topPatchID, bottomPatchID)
312  // - per edge, per face on edge the side patch (edgePatches)
313  labelListList edgePatches(extrudedPatch_.nEdges());
314  forAll(edgePatches, edgeI)
315  {
316  const labelList& eFaces = extrudedPatch_.edgeFaces()[edgeI];
317 
318  if (eFaces.size() != 2 || nonManifoldEdge[edgeI])
319  {
320  edgePatches[edgeI].setSize(eFaces.size(), sidePatchID);
321  }
322  }
323 
324  polyTopoChange meshMod(regionPatches.size());
325 
326  extruder.setRefinement
327  (
328  firstDisp, // first displacement
329  model_().expansionRatio(),
330  model_().nLayers(), // nLayers
331  labelList(extrudedPatch_.size(), topPatchID),
332  labelList(extrudedPatch_.size(), bottomPatchID),
333  edgePatches,
334  meshMod
335  );
336 
337  autoPtr<mapPolyMesh> map = meshMod.changeMesh
338  (
339  *this, // mesh to change
340  false // inflate
341  );
342 
343  // Update numbering on extruder.
344  extruder.updateMesh(map);
345 
346  this->setInstance(this->thisDb().time().constant());
347  this->write();
348  }
349 }
350 
351 
352 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
353 
355 {}
356 
357 
358 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
359 
360 } // End namespace Foam
361 
362 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
void clearOut()
Clear all geometry and addressing.
Definition: fvMesh.C:240
static void calcPointRegions(const globalMeshData &globalData, const primitiveFacePatch &patch, const PackedBoolList &nonManifoldEdge, const bool syncNonCollocated, faceList &pointGlobalRegions, faceList &pointLocalRegions, labelList &localToGlobalRegion)
Helper: calculate point regions. The point region is the.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
void removeFvBoundary()
Remove boundary patches. Warning: fvPatchFields hold ref to.
Definition: fvMesh.C:473
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
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:808
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
extrudePatchMesh(const fvMesh &, const fvPatch &, const dictionary &, const word)
Construct from mesh, patch and dictionary.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
void addFvPatches(const List< polyPatch *> &, const bool validBoundary=true)
Add boundary patches. Constructor helper.
Definition: fvMesh.C:455
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
wordList patchTypes(nPatches)
List< face > faceList
Definition: faceListFwd.H:43
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
label nInternalEdges() const
Number of internal edges.
bool readBool(Istream &)
Definition: boolIO.C:60
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
virtual bool write(const bool write=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:1093
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
virtual const objectRegistry & thisDb() const
Return the object registry - resolve conflict polyMesh/lduMesh.
Definition: fvMesh.H:245
static void mapCombineScatter(const List< commsStruct > &comms, Container &Values, const int tag, const label comm)
Scatter data. Reverse of combineGather.
A list of faces which address into the list of points.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
dynamicFvMesh & mesh
static autoPtr< extrudeModel > New(const dictionary &)
Select null constructed.
A class for handling words, derived from string.
Definition: word.H:59
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:138
wordList patchNames(nPatches)
Creates mesh by extruding a patch.
const Field< PointType > & points() const
Return reference to global points.
const labelListList & edgeFaces() const
Return edge-face addressing.
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1394
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const zero Zero
Definition: zero.H:97
static void mapCombineGather(const List< commsStruct > &comms, Container &Values, const CombineOp &cop, const int tag, const label comm)
const Field< PointType > & faceNormals() const
Return face normals for patch.
label nEdges() const
Return number of edges in patch.
void updateMesh(const mapPolyMesh &)
Update any locally stored mesh information.
defineTypeNameAndDebug(combustionModel, 0)
labelList f(nPoints)
void setInstance(const fileName &)
Set the instance for mesh files.
Definition: polyMeshIO.C:32
void setSize(const label)
Reset size of List.
Definition: List.C:281
A bit-packed bool list.
label patchi
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
const List< FaceType > & localFaces() const
Return patch faces addressing into local point list.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Direct mesh changes based on v1.3 polyTopoChange syntax.
void set(entry *)
Assign a new entry, overwrite any existing entry.
Definition: dictionary.C:1271
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
void setRefinement(const pointField &firstLayerThickness, const scalar expansionRatio, const label nLayers, const labelList &topPatchID, const labelList &bottomPatchID, const labelListList &extrudeEdgePatches, polyTopoChange &meshMod)
Play commands into polyTopoChange to create layer mesh.
virtual ~extrudePatchMesh()
Destructor.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844
A HashTable to objects of type <T> with a label key.
Definition: Map.H:49