channelIndex.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 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 "channelIndex.H"
27 #include "boolList.H"
28 #include "syncTools.H"
29 #include "OFstream.H"
30 #include "meshTools.H"
31 #include "Time.H"
32 #include "SortableList.H"
33 
34 // * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  template<>
39  const char* Foam::NamedEnum
40  <
42  3
43  >::names[] =
44  {
45  "x",
46  "y",
47  "z"
48  };
49 }
50 
52  Foam::channelIndex::vectorComponentsNames_;
53 
54 
55 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
56 
57 // Determines face blocking
58 void Foam::channelIndex::walkOppositeFaces
59 (
60  const polyMesh& mesh,
61  const labelList& startFaces,
62  boolList& blockedFace
63 )
64 {
65  const cellList& cells = mesh.cells();
66  const faceList& faces = mesh.faces();
67  label nBnd = mesh.nFaces() - mesh.nInternalFaces();
68 
69  DynamicList<label> frontFaces(startFaces);
70  forAll(frontFaces, i)
71  {
72  label facei = frontFaces[i];
73  blockedFace[facei] = true;
74  }
75 
76  while (returnReduce(frontFaces.size(), sumOp<label>()) > 0)
77  {
78  // Transfer across.
79  boolList isFrontBndFace(nBnd, false);
80  forAll(frontFaces, i)
81  {
82  label facei = frontFaces[i];
83 
84  if (!mesh.isInternalFace(facei))
85  {
86  isFrontBndFace[facei-mesh.nInternalFaces()] = true;
87  }
88  }
89  syncTools::swapBoundaryFaceList(mesh, isFrontBndFace);
90 
91  // Add
92  forAll(isFrontBndFace, i)
93  {
94  label facei = mesh.nInternalFaces()+i;
95  if (isFrontBndFace[i] && !blockedFace[facei])
96  {
97  blockedFace[facei] = true;
98  frontFaces.append(facei);
99  }
100  }
101 
102  // Transfer across cells
103  DynamicList<label> newFrontFaces(frontFaces.size());
104 
105  forAll(frontFaces, i)
106  {
107  label facei = frontFaces[i];
108 
109  {
110  const cell& ownCell = cells[mesh.faceOwner()[facei]];
111 
112  label oppositeFacei = ownCell.opposingFaceLabel(facei, faces);
113 
114  if (oppositeFacei == -1)
115  {
117  << "Face:" << facei << " owner cell:" << ownCell
118  << " is not a hex?" << abort(FatalError);
119  }
120  else
121  {
122  if (!blockedFace[oppositeFacei])
123  {
124  blockedFace[oppositeFacei] = true;
125  newFrontFaces.append(oppositeFacei);
126  }
127  }
128  }
129 
130  if (mesh.isInternalFace(facei))
131  {
132  const cell& neiCell = mesh.cells()[mesh.faceNeighbour()[facei]];
133 
134  label oppositeFacei = neiCell.opposingFaceLabel(facei, faces);
135 
136  if (oppositeFacei == -1)
137  {
139  << "Face:" << facei << " neighbour cell:" << neiCell
140  << " is not a hex?" << abort(FatalError);
141  }
142  else
143  {
144  if (!blockedFace[oppositeFacei])
145  {
146  blockedFace[oppositeFacei] = true;
147  newFrontFaces.append(oppositeFacei);
148  }
149  }
150  }
151  }
152 
153  frontFaces.transfer(newFrontFaces);
154  }
155 }
156 
157 
158 // Calculate regions.
159 void Foam::channelIndex::calcLayeredRegions
160 (
161  const polyMesh& mesh,
162  const labelList& startFaces
163 )
164 {
165  boolList blockedFace(mesh.nFaces(), false);
166  walkOppositeFaces
167  (
168  mesh,
169  startFaces,
170  blockedFace
171  );
172 
173 
174  if (false)
175  {
176  OFstream str(mesh.time().path()/"blockedFaces.obj");
177  label vertI = 0;
178  forAll(blockedFace, facei)
179  {
180  if (blockedFace[facei])
181  {
182  const face& f = mesh.faces()[facei];
183  forAll(f, fp)
184  {
185  meshTools::writeOBJ(str, mesh.points()[f[fp]]);
186  }
187  str<< 'f';
188  forAll(f, fp)
189  {
190  str << ' ' << vertI+fp+1;
191  }
192  str << nl;
193  vertI += f.size();
194  }
195  }
196  }
197 
198 
199  // Do analysis for connected regions
200  cellRegion_.reset(new regionSplit(mesh, blockedFace));
201 
202  Info<< "Detected " << cellRegion_().nRegions() << " layers." << nl << endl;
203 
204  // Sum number of entries per region
205  regionCount_ = regionSum(scalarField(mesh.nCells(), 1.0));
206 
207  // Average cell centres to determine ordering.
208  pointField regionCc
209  (
210  regionSum(mesh.cellCentres())
211  / regionCount_
212  );
213 
214  SortableList<scalar> sortComponent(regionCc.component(dir_));
215 
216  sortMap_ = sortComponent.indices();
217 
218  y_ = sortComponent;
219 
220  if (symmetric_)
221  {
222  y_.setSize(cellRegion_().nRegions()/2);
223  }
224 }
225 
226 
227 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
228 
229 Foam::channelIndex::channelIndex
230 (
231  const polyMesh& mesh,
232  const dictionary& dict
233 )
234 :
235  symmetric_(readBool(dict.lookup("symmetric"))),
236  dir_(vectorComponentsNames_.read(dict.lookup("component")))
237 {
238  const polyBoundaryMesh& patches = mesh.boundaryMesh();
239 
240  const wordList patchNames(dict.lookup("patches"));
241 
242  label nFaces = 0;
243 
244  forAll(patchNames, i)
245  {
246  const label patchi = patches.findPatchID(patchNames[i]);
247 
248  if (patchi == -1)
249  {
251  << "Illegal patch " << patchNames[i]
252  << ". Valid patches are " << patches.name()
253  << exit(FatalError);
254  }
255 
256  nFaces += patches[patchi].size();
257  }
258 
259  labelList startFaces(nFaces);
260  nFaces = 0;
261 
262  forAll(patchNames, i)
263  {
264  const polyPatch& pp = patches[patchNames[i]];
265 
266  forAll(pp, j)
267  {
268  startFaces[nFaces++] = pp.start()+j;
269  }
270  }
271 
272  // Calculate regions.
273  calcLayeredRegions(mesh, startFaces);
274 }
275 
276 
277 Foam::channelIndex::channelIndex
278 (
279  const polyMesh& mesh,
280  const labelList& startFaces,
281  const bool symmetric,
282  const direction dir
283 )
284 :
285  symmetric_(symmetric),
286  dir_(dir)
287 {
288  // Calculate regions.
289  calcLayeredRegions(mesh, startFaces);
290 }
291 
292 
293 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
uint8_t direction
Definition: direction.H:46
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
List< face > faceList
Definition: faceListFwd.H:43
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
bool readBool(Istream &)
Definition: boolIO.C:60
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:52
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
bool read(const char *, int32_t &)
Definition: int32IO.C:85
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
stressControl lookup("compactNormalStress") >> compactNormalStress
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
wordList patchNames(nPatches)
List< label > labelList
A List of labels.
Definition: labelList.H:56
errorManip< error > abort(error &err)
Definition: errorManip.H:131
static const char nl
Definition: Ostream.H:262
List< word > wordList
A List of words.
Definition: fileName.H:54
label patchi
messageStream Info
components
Component labeling enumeration.
Definition: Vector.H:75
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &l)
Swap coupled boundary face values.
Definition: syncTools.H:430
List< cell > cellList
list of cells
Definition: cellList.H:42
Namespace for OpenFOAM.