blockMeshCheck.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 "blockMesh.H"
27 
28 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
29 
30 void Foam::blockMesh::checkBlockMesh(const polyMesh& bm) const
31 {
32  if (verboseOutput)
33  {
34  Info<< nl << "Check topology" << endl;
35  }
36 
37  bool ok = true;
38 
39  const pointField& points = bm.points();
40  const faceList& faces = bm.faces();
41  const cellList& cells = bm.cells();
42  const polyPatchList& patches = bm.boundaryMesh();
43 
44  label nBoundaryFaces = 0;
45  forAll(cells, celli)
46  {
47  nBoundaryFaces += cells[celli].nFaces();
48  }
49 
50  nBoundaryFaces -= 2*bm.nInternalFaces();
51 
52  label nDefinedBoundaryFaces = 0;
53  forAll(patches, patchi)
54  {
55  nDefinedBoundaryFaces += patches[patchi].size();
56  }
57 
58 
59  if (verboseOutput)
60  {
61  Info<< nl << tab << "Basic statistics" << nl
62  << tab << tab << "Number of internal faces : "
63  << bm.nInternalFaces() << nl
64  << tab << tab << "Number of boundary faces : "
65  << nBoundaryFaces << nl
66  << tab << tab << "Number of defined boundary faces : "
67  << nDefinedBoundaryFaces << nl
68  << tab << tab << "Number of undefined boundary faces : "
69  << nBoundaryFaces - nDefinedBoundaryFaces << nl;
70 
71  if ((nBoundaryFaces - nDefinedBoundaryFaces) > 0)
72  {
73  Info<< tab << tab << tab
74  << "(Warning : only leave undefined the front and back planes "
75  << "of 2D planar geometries!)" << endl;
76  }
77 
78  Info<< tab << "Checking patch -> block consistency" << endl;
79  }
80 
81 
82  forAll(patches, patchi)
83  {
84  const faceList& Patch = patches[patchi];
85 
86  forAll(Patch, patchFacei)
87  {
88  const face& patchFace = Patch[patchFacei];
89  bool patchFaceOK = false;
90 
91  forAll(cells, celli)
92  {
93  const labelList& cellFaces = cells[celli];
94 
95  forAll(cellFaces, cellFacei)
96  {
97  if (patchFace == faces[cellFaces[cellFacei]])
98  {
99  patchFaceOK = true;
100 
101  if
102  (
103  (
104  patchFace.normal(points)
105  & faces[cellFaces[cellFacei]].normal(points)
106  ) < 0.0
107  )
108  {
109  Info<< tab << tab
110  << "Face " << patchFacei
111  << " of patch " << patchi
112  << " (" << patches[patchi].name() << ")"
113  << " points inwards"
114  << endl;
115 
116  ok = false;
117  }
118  }
119  }
120  }
121 
122  if (!patchFaceOK)
123  {
124  Info<< tab << tab
125  << "Face " << patchFacei
126  << " of patch " << patchi
127  << " (" << patches[patchi].name() << ")"
128  << " does not match any block faces" << endl;
129 
130  ok = false;
131  }
132  }
133  }
134 
135  if (verboseOutput)
136  {
137  Info<< endl;
138  }
139 
140  if (!ok)
141  {
143  << "Block mesh topology incorrect, stopping mesh generation!"
144  << exit(FatalError);
145  }
146 }
147 
148 
149 bool Foam::blockMesh::blockLabelsOK
150 (
151  const label blockLabel,
152  const pointField& points,
153  const cellShape& blockShape
154 ) const
155 {
156  bool ok = true;
157 
158  forAll(blockShape, blockI)
159  {
160  if (blockShape[blockI] < 0)
161  {
162  ok = false;
163 
165  << "out-of-range point label " << blockShape[blockI]
166  << " (min = 0"
167  << ") in block " << blockLabel << endl;
168  }
169  else if (blockShape[blockI] >= points.size())
170  {
171  ok = false;
172 
174  << "out-of-range point label " << blockShape[blockI]
175  << " (max = " << points.size() - 1
176  << ") in block " << blockLabel << endl;
177  }
178  }
179 
180  return ok;
181 }
182 
183 
184 bool Foam::blockMesh::patchLabelsOK
185 (
186  const label patchLabel,
187  const pointField& points,
188  const faceList& patchFaces
189 ) const
190 {
191  bool ok = true;
192 
193  forAll(patchFaces, facei)
194  {
195  const labelList& f = patchFaces[facei];
196 
197  forAll(f, fp)
198  {
199  if (f[fp] < 0)
200  {
201  ok = false;
202 
204  << "out-of-range point label " << f[fp]
205  << " (min = 0"
206  << ") on patch " << patchLabel
207  << ", face " << facei << endl;
208  }
209  else if (f[fp] >= points.size())
210  {
211  ok = false;
212 
214  << "out-of-range point label " << f[fp]
215  << " (max = " << points.size() - 1
216  << ") on patch " << patchLabel
217  << ", face " << facei << endl;
218 
219  }
220  }
221  }
222 
223  return ok;
224 }
225 
226 
227 // ************************************************************************* //
PtrList< polyPatch > polyPatchList
container classes for polyPatch
Definition: polyPatchList.H:45
#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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const char tab
Definition: Ostream.H:261
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
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const char nl
Definition: Ostream.H:262
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
messageStream Info
List< cell > cellList
list of cells
Definition: cellList.H:42