mappedExtrudedPatchBase.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) 2022-2023 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 
27 #include "LayerInfoData.H"
28 #include "PointEdgeLayerInfoData.H"
29 #include "FaceCellWave.H"
30 #include "PointEdgeWave.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
37 
38  template<>
39  inline bool contiguous<LayerInfoData<Pair<vector>>>()
40  {
41  return true;
42  }
43 }
44 
45 
46 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
47 
50 {
51  if (isExtrudedRegion_)
52  {
53  if (!bottomFaceAreasPtr_.valid())
54  {
55  const polyMesh& mesh = patch_.boundaryMesh().mesh();
56  const polyPatch& pp = patch_;
57 
58  // If this is the extruded region we need to work out what the
59  // corresponding areas and centres are on the bottom patch. We do
60  // this by waving these values across the layers.
61 
62  // Initialise layer data on the patch faces
63  labelList initialFaces1(pp.size());
64  List<layerInfo> initialFaceInfo1(pp.size());
65  forAll(pp, ppFacei)
66  {
67  initialFaces1[ppFacei] = pp.start() + ppFacei;
68  initialFaceInfo1[ppFacei] = layerInfo(0, -1);
69  }
70 
71  // Wave across the mesh layers
72  List<layerInfo> faceInfo1(mesh.nFaces());
73  List<layerInfo> cellInfo1(mesh.nCells());
74  FaceCellWave<layerInfo> wave1(mesh, faceInfo1, cellInfo1);
75  wave1.setFaceInfo(initialFaces1, initialFaceInfo1);
76  const label nIterations1 =
77  wave1.iterate(mesh.globalData().nTotalCells() + 1);
78 
79  // Count how many opposite faces the wave ended on
80  label nInitialFaces2 = 0;
81  for
82  (
83  label facei = mesh.nInternalFaces();
84  facei < mesh.nFaces();
85  ++ facei
86  )
87  {
88  if
89  (
90  faceInfo1[facei].valid(wave1.data())
91  && faceInfo1[facei].faceLayer() == nIterations1
92  )
93  {
94  nInitialFaces2 ++;
95  }
96  }
97 
98  // Initialise data on the opposite faces. Store the area and centre.
99  labelList initialFaces2(nInitialFaces2);
100  List<LayerInfoData<Pair<vector>>> initialFaceInfo2(nInitialFaces2);
101  label initialFace2i = 0;
102  for
103  (
104  label facei = mesh.nInternalFaces();
105  facei < mesh.nFaces();
106  ++ facei
107  )
108  {
109  if
110  (
111  faceInfo1[facei].valid(wave1.data())
112  && faceInfo1[facei].faceLayer() == nIterations1
113  )
114  {
115  initialFaces2[initialFace2i] = facei;
116  initialFaceInfo2[initialFace2i] =
118  (
119  0,
120  -1,
122  (
123  mesh.faceAreas()[facei],
124  mesh.faceCentres()[facei]
125  )
126  );
127  initialFace2i ++;
128  }
129  }
130 
131  // Wave back across the mesh layers
132  List<LayerInfoData<Pair<vector>>> faceInfo2(mesh.nFaces());
133  List<LayerInfoData<Pair<vector>>> cellInfo2(mesh.nCells());
135  (
136  mesh,
137  initialFaces2,
138  initialFaceInfo2,
139  faceInfo2,
140  cellInfo2,
141  mesh.globalData().nTotalCells() + 1
142  );
143 
144  // Unpack into this patch's bottom face areas and centres. Note
145  // that the face area needs flipping as it relates to a patch on
146  // the other side of the extruded region.
147  bottomFaceAreasPtr_.set(new vectorField(pp.size()));
148  bottomFaceCentresPtr_.set(new pointField(pp.size()));
149  forAll(pp, ppFacei)
150  {
151  const LayerInfoData<Pair<vector>>& info =
152  faceInfo2[pp.start() + ppFacei];
153 
154  static nil td;
155 
156  if (!info.valid(td))
157  {
159  << "Mesh \"" << mesh.name()
160  << "\" is not layered from the extruded patch "
161  << "\"" << pp.name() << "\"" << exit(FatalError);
162  }
163 
164  bottomFaceAreasPtr_()[ppFacei] = - info.data().first();
165  bottomFaceCentresPtr_()[ppFacei] = info.data().second();
166  }
167  }
168 
169  return bottomFaceAreasPtr_();
170  }
171 
173 }
174 
175 
178 {
179  if (isExtrudedRegion_)
180  {
181  if (!bottomFaceCentresPtr_.valid())
182  {
183  patchFaceAreas();
184  }
185 
186  return bottomFaceCentresPtr_();
187  }
188 
190 }
191 
192 
195 {
196  if (isExtrudedRegion_)
197  {
198  if (!bottomLocalPointsPtr_.valid())
199  {
200  const polyMesh& mesh = patch_.boundaryMesh().mesh();
201  const polyPatch& pp = patch_;
202 
203  // If this is the extruded region we need to work out what the
204  // corresponding points are on the bottom patch. We do this by
205  // waving these location across the layers.
206 
207  // Initialise layer data on the patch points
208  labelList initialPoints1(pp.nPoints());
209  List<pointEdgeLayerInfo> initialPointInfo1(pp.nPoints());
210  forAll(pp.meshPoints(), ppPointi)
211  {
212  initialPoints1[ppPointi] = pp.meshPoints()[ppPointi];
213  initialPointInfo1[ppPointi] = pointEdgeLayerInfo(0);
214  }
215 
216  // Wave across the mesh layers
217  List<pointEdgeLayerInfo> pointInfo1(mesh.nPoints());
218  List<pointEdgeLayerInfo> edgeInfo1(mesh.nEdges());
220  (
221  mesh,
222  pointInfo1,
223  edgeInfo1
224  );
225  wave1.setPointInfo(initialPoints1, initialPointInfo1);
226  const label nIterations1 =
227  wave1.iterate(mesh.globalData().nTotalPoints() + 1);
228 
229  if (debug)
230  {
231  pointScalarField pointLayer
232  (
234  (
235  typedName("pointLayer"),
236  pointMesh::New(mesh),
238  )
239  );
240  forAll(pointInfo1, pointi)
241  {
242  pointLayer[pointi] =
243  pointInfo1[pointi].valid(wave1.data())
244  ? pointInfo1[pointi].pointLayer()
245  : -1;
246  }
247  pointLayer.write();
248  }
249 
250  // Count how many opposite points the wave ended on
251  label nInitialPoints2 = 0;
252  forAll(pointInfo1, pointi)
253  {
254  if
255  (
256  pointInfo1[pointi].valid(wave1.data())
257  && pointInfo1[pointi].pointLayer() == nIterations1
258  )
259  {
260  nInitialPoints2 ++;
261  }
262  }
263 
264  // Initialise data on the opposite points. Store the position.
265  labelList initialPoints2(nInitialPoints2);
267  initialPointInfo2(nInitialPoints2);
268  label initialPoint2i = 0;
269  forAll(pointInfo1, pointi)
270  {
271  if
272  (
273  pointInfo1[pointi].valid(wave1.data())
274  && pointInfo1[pointi].pointLayer() == nIterations1
275  )
276  {
277  initialPoints2[initialPoint2i] = pointi;
278  initialPointInfo2[initialPoint2i] =
279  PointEdgeLayerInfoData<point>(0, mesh.points()[pointi]);
280  initialPoint2i ++;
281  }
282  }
283 
284  // Wave back across the mesh layers
285  List<PointEdgeLayerInfoData<point>> pointInfo2(mesh.nPoints());
286  List<PointEdgeLayerInfoData<point>> edgeInfo2(mesh.nEdges());
288  (
289  mesh,
290  initialPoints2,
291  initialPointInfo2,
292  pointInfo2,
293  edgeInfo2,
294  mesh.globalData().nTotalCells() + 1
295  );
296 
297  // Unpack into this patch's bottom local points
298  bottomLocalPointsPtr_.set(new pointField(pp.nPoints()));
299  forAll(pp.meshPoints(), ppPointi)
300  {
301  const PointEdgeLayerInfoData<point>& info =
302  pointInfo2[pp.meshPoints()[ppPointi]];
303 
304  static nil td;
305 
306  if (!info.valid(td))
307  {
309  << "Mesh \"" << mesh.name()
310  << "\" is not layered from the extruded patch "
311  << "\"" << pp.name() << "\"" << exit(FatalError);
312  }
313 
314  bottomLocalPointsPtr_()[ppPointi] = info.data();
315  }
316 
317  if (debug)
318  {
319  pointVectorField pointOffset
320  (
322  (
323  typedName("pointOffset"),
324  pointMesh::New(mesh),
326  )
327  );
328  forAll(pp.meshPoints(), ppPointi)
329  {
330  pointOffset[pp.meshPoints()[ppPointi]] =
331  bottomLocalPointsPtr_()[ppPointi]
332  - pp.localPoints()[ppPointi];
333  }
334  pointOffset.write();
335  }
336  }
337 
338  return bottomLocalPointsPtr_();
339  }
340 
342 }
343 
344 
345 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
346 
348 :
349  mappedPatchBase(pp),
350  isExtrudedRegion_(false)
351 {}
352 
353 
355 (
356  const polyPatch& pp,
357  const word& nbrRegionName,
358  const word& nbrPatchName,
359  const bool isExtrudedRegion,
361 )
362 :
363  mappedPatchBase(pp, nbrRegionName, nbrPatchName, transform),
364  isExtrudedRegion_(isExtrudedRegion)
365 {}
366 
367 
369 (
370  const polyPatch& pp,
371  const dictionary& dict,
372  const bool defaultTransformIsNone
373 )
374 :
375  mappedPatchBase(pp, dict, defaultTransformIsNone),
376  isExtrudedRegion_(dict.lookup<bool>("isExtrudedRegion"))
377 {}
378 
379 
381 (
382  const polyPatch& pp,
383  const mappedExtrudedPatchBase& mepb
384 )
385 :
386  mappedPatchBase(pp, mepb),
387  isExtrudedRegion_(mepb.isExtrudedRegion_)
388 {}
389 
390 
391 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
392 
394 {}
395 
396 
397 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
398 
400 {
402  if (reMapAfterMove_)
403  {
404  bottomFaceCentresPtr_.clear();
405  }
406 }
407 
408 
410 {
412  writeEntry(os, "isExtrudedRegion", isExtrudedRegion_);
413 }
414 
415 
416 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
static pointMesh & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
Wave propagation of information through grid. Every iteration information goes through one layer of c...
Definition: FaceCellWave.H:79
const TrackingData & data() const
Additional data to be passed into container.
Definition: FaceCellWave.H:308
void setFaceInfo(const labelList &changedFaces, const List< Type > &changedFacesInfo)
Set initial changed faces.
Definition: FaceCellWave.C:260
virtual label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
Definition: FaceCellWave.C:951
static tmp< GeometricField< Type, pointPatchField, pointMesh > > New(const word &name, const Internal &, const PtrList< pointPatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
const word & name() const
Return name.
Definition: IOobject.H:310
Class to be used with FaceCellWave which enumerates layers of cells and transports data through those...
Definition: LayerInfoData.H:66
const Type & data() const
Return the data.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: Pair.H:65
Class to be used with PointEdgeWave which enumerates layers of points.
const Type & data() const
Return the data.
Wave propagation of information through grid. Every iteration information goes through one layer of e...
Definition: PointEdgeWave.H:89
const TrackingData & data() const
Additional data to be passed into container.
void setPointInfo(const labelList &changedPoints, const List< Type > &changedPointsInfo)
Copy initial data into allPointInfo_.
label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
label nPoints() const
Return number of points supporting patch faces.
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
const Field< PointType > & localPoints() const
Return pointField of points in patch.
label size() const
Return the number of elements in the UList.
Definition: UListI.H:311
Cyclic plane transformation.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
label nTotalPoints() const
Return total number of points in decomposed mesh. Not.
label nTotalCells() const
Return total number of cells in decomposed mesh.
Class to be used with FaceCellWave which enumerates layers of cells.
Definition: layerInfo.H:60
bool valid(TrackingData &td) const
Check whether the layerInfo has been changed at all or still.
Definition: layerInfoI.H:88
Engine which provides mapping between two patches and which offsets the geometry to account for the e...
mappedExtrudedPatchBase(const polyPatch &)
Construct from patch.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual tmp< vectorField > patchFaceAreas() const
Get the face-areas for the patch.
virtual tmp< pointField > patchLocalPoints() const
Get the local points for the patch.
virtual ~mappedExtrudedPatchBase()
Destructor.
virtual tmp< pointField > patchFaceCentres() const
Get the face-centres for the patch.
virtual void clearOut()
Clear out data on mesh change.
const polyPatch & patch_
Patch to map to.
Engine and base class for poly patches which provides interpolative mapping between two globally conf...
virtual void write(Ostream &) const
Write as a dictionary.
virtual tmp< vectorField > patchFaceAreas() const
Get the face-areas for this patch.
virtual tmp< pointField > patchLocalPoints() const
Get the local points for this patch.
virtual tmp< pointField > patchFaceCentres() const
Get the face-centres for this patch.
virtual void clearOut()
Clear out data on mesh change.
A zero-sized class without any storage. Used, for example, in HashSet.
Definition: nil.H:59
const word & name() const
Return name.
Class to be used with PointEdgeWave which enumerates layers of points.
bool valid(TrackingData &td) const
Check whether info has been changed at all or.
const polyMesh & mesh() const
Return the mesh reference.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1515
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:404
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1313
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:280
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:270
label nEdges() const
label nCells() const
label nPoints() const
const vectorField & faceCentres() const
label nInternalFaces() const
label nFaces() const
const vectorField & faceAreas() const
virtual bool write(const bool write=true) const
Write using setting from DB.
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
bool valid(const PtrList< ModelType > &l)
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const zero Zero
Definition: zero.H:97
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 dimensionSet dimless
const dimensionSet dimLength
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
defineTypeNameAndDebug(combustionModel, 0)
word typedName(Name name)
Return the name of the object within the given type.
Definition: typeInfo.H:176
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Field< vector > vectorField
Specialisation of Field<T> for vector.
error FatalError
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dictionary dict