wallLayerCells.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-2018 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 "wallLayerCells.H"
27 #include "DynamicList.H"
28 #include "MeshWave.H"
29 #include "wallNormalInfo.H"
30 #include "OFstream.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(wallLayerCells, 0);
37 }
38 
39 
40 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
41 
42 bool Foam::wallLayerCells::usesCoupledPatch(const label celli) const
43 {
44  const polyBoundaryMesh& patches = mesh().boundaryMesh();
45 
46  const cell& cFaces = mesh().cells()[celli];
47 
48  forAll(cFaces, cFacei)
49  {
50  label facei = cFaces[cFacei];
51 
52  label patchID = patches.whichPatch(facei);
53 
54  if ((patchID >= 0) && (patches[patchID].coupled()))
55  {
56  return true;
57  }
58  }
59  return false;
60 }
61 
62 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
63 
65 (
66  const polyMesh& mesh,
67  const List<word>& patchNames,
68  const label nLayers
69 )
70 :
71  edgeVertex(mesh),
73 {
74  // Find out cells connected to walls.
75 
76  const polyPatchList& patches = mesh.boundaryMesh();
77 
78  // Make map from name to local patch ID
79  HashTable<label> patchNameToIndex(patches.size());
80 
81  forAll(patches, patchi)
82  {
83  patchNameToIndex.insert(patches[patchi].name(), patchi);
84  }
85 
86 
87  // Count size of walls to set
88  label nWalls = 0;
89 
90  forAll(patchNames, patchNameI)
91  {
92  const word& name = patchNames[patchNameI];
93 
94  if (patchNameToIndex.found(name))
95  {
96  label patchi = patchNameToIndex[name];
97 
98  nWalls += patches[patchi].size();
99  }
100  }
101 
102  // Allocate storage for start of wave on faces
103  List<wallNormalInfo> changedFacesInfo(nWalls);
104  labelList changedFaces(nWalls);
105 
106  // Fill changedFaces info
107  label nChangedFaces = 0;
108 
109  forAll(patchNames, patchNameI)
110  {
111  const word& name = patchNames[patchNameI];
112 
113  if (patchNameToIndex.found(name))
114  {
115  label patchi = patchNameToIndex[name];
116 
117  const polyPatch& pp = patches[patchi];
118 
119  forAll(pp, patchFacei)
120  {
121  label meshFacei = pp.start() + patchFacei;
122 
123  changedFaces[nChangedFaces] = meshFacei;
124 
125  // Set transported information to the wall normal.
126  const vector& norm = pp.faceNormals()[patchFacei];
127 
128  changedFacesInfo[nChangedFaces] = wallNormalInfo(norm);
129 
130  nChangedFaces++;
131  }
132  }
133  }
134 
135 
136  // Do a wave of nLayers, transporting the index in patchNames
137  // (cannot use local patchIDs since we might get info from neighbouring
138  // processor)
139 
140  MeshWave<wallNormalInfo> regionCalc
141  (
142  mesh,
143  changedFaces,
144  changedFacesInfo,
145  0
146  );
147 
148  regionCalc.iterate(nLayers);
149 
150 
151  // Now regionCalc should hold info on faces that are reachable from
152  // changedFaces within nLayers iterations. We use face info since that is
153  // guaranteed to be consistent across processor boundaries.
154 
155  const List<wallNormalInfo>& faceInfo = regionCalc.allFaceInfo();
156 
157  if (debug)
158  {
159  InfoInFunction << "Dumping selected faces to selectedFaces.obj" << endl;
160 
161  OFstream fcStream("selectedFaces.obj");
162 
163  label vertI = 0;
164 
165  forAll(faceInfo, facei)
166  {
167  const wallNormalInfo& info = faceInfo[facei];
168 
169  if (info.valid(regionCalc.data()))
170  {
171  const face& f = mesh.faces()[facei];
172 
173  point mid(0.0, 0.0, 0.0);
174 
175  forAll(f, fp)
176  {
177  mid += mesh.points()[f[fp]];
178  }
179  mid /= f.size();
180 
181  fcStream
182  << "v " << mid.x() << ' ' << mid.y() << ' ' << mid.z()
183  << endl;
184  vertI++;
185 
186  point end(mid + info.normal());
187 
188  fcStream
189  << "v " << end.x() << ' ' << end.y() << ' ' << end.z()
190  << endl;
191  vertI++;
192 
193  fcStream << "l " << vertI << ' ' <<vertI-1 << endl;
194  }
195  }
196  }
197 
198 
199  //
200  // Copy meshWave information to List<refineCell>
201  //
202 
203  // Estimate size
204 
205  DynamicList<refineCell> refineCells(3*nWalls);
206 
207  const List<wallNormalInfo>& cellInfo = regionCalc.allCellInfo();
208 
209  forAll(cellInfo, celli)
210  {
211  const wallNormalInfo& info = cellInfo[celli];
212 
213  if (info.valid(regionCalc.data()) && !usesCoupledPatch(celli))
214  {
215  refineCells.append(refineCell(celli, info.normal()));
216  }
217  }
218 
219  // Transfer refineCells storage to this.
220  transfer(refineCells);
221 }
222 
223 
224 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
wallLayerCells(const polyMesh &mesh, const List< word > &patchNames, const label nLayers)
Construct from components.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
Holds information regarding nearest wall point. Used in wall refinement.
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
Output to file stream.
Definition: OFstream.H:82
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const vector & normal() const
const Cmpt & z() const
Definition: VectorI.H:87
const cellList & cells() const
Holds information regarding type of cell. Used in inside/outside determination in cellClassification...
Definition: cellInfo.H:63
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
const TrackingData & data() const
Additional data to be passed into container.
Definition: MeshWave.H:128
const Cmpt & y() const
Definition: VectorI.H:81
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
const List< Type > & allFaceInfo() const
Get allFaceInfo.
Definition: MeshWave.H:116
dynamicFvMesh & mesh
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
A class for handling words, derived from string.
Definition: word.H:59
FaceCellWave plus data.
Definition: MeshWave.H:56
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
Combines edge or vertex in single label. Used to specify cuts across cell circumference.
Definition: edgeVertex.H:52
label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
Definition: MeshWave.H:135
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
const Field< PointType > & faceNormals() const
Return face normals for patch.
Container with cells to refine. Refinement given as single direction.
Definition: refineCell.H:56
const Cmpt & x() const
Definition: VectorI.H:75
defineTypeNameAndDebug(combustionModel, 0)
labelList f(nPoints)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const List< Type > & allCellInfo() const
Get allCellInfo.
Definition: MeshWave.H:122
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
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
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:309
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Namespace for OpenFOAM.
#define InfoInFunction
Report an information message using Foam::Info.