polyMeshInitMesh.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 "polyMesh.H"
27 
28 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
29 
30 void Foam::polyMesh::initMesh()
31 {
32  if (debug)
33  {
34  InfoInFunction << "initialising primitiveMesh" << endl;
35  }
36 
37  // For backward compatibility check if the neighbour array is the same
38  // length as the owner and shrink to remove the -1s padding
39  if (neighbour_.size() == owner_.size())
40  {
42 
43  forAll(neighbour_, facei)
44  {
45  if (neighbour_[facei] == -1)
46  {
47  break;
48  }
49  else
50  {
52  }
53  }
54 
55  neighbour_.setSize(nInternalFaces);
56  }
57 
58  label nCells = -1;
59 
60  forAll(owner_, facei)
61  {
62  if (owner_[facei] < 0)
63  {
65  << "Illegal cell label " << owner_[facei]
66  << " in neighbour addressing for face " << facei
67  << exit(FatalError);
68  }
69  nCells = max(nCells, owner_[facei]);
70  }
71 
72  // The neighbour array may or may not be the same length as the owner
73  forAll(neighbour_, facei)
74  {
75  if (neighbour_[facei] < 0)
76  {
78  << "Illegal cell label " << neighbour_[facei]
79  << " in neighbour addressing for face " << facei
80  << exit(FatalError);
81  }
82  nCells = max(nCells, neighbour_[facei]);
83  }
84 
85  nCells++;
86 
87  // Reset the primitiveMesh with the sizes of the primitive arrays
89  (
90  points_.size(),
91  neighbour_.size(),
92  owner_.size(),
93  nCells
94  );
95 
96  string meshInfo =
97  "nPoints:" + Foam::name(nPoints())
98  + " nCells:" + Foam::name(this->nCells())
99  + " nFaces:" + Foam::name(nFaces())
100  + " nInternalFaces:" + Foam::name(nInternalFaces());
101 
102  owner_.note() = meshInfo;
103  neighbour_.note() = meshInfo;
104 }
105 
106 
107 void Foam::polyMesh::initMesh(cellList& c)
108 {
109  if (debug)
110  {
111  InfoInFunction << "Calculating owner-neighbour arrays" << endl;
112  }
113 
114  owner_.setSize(faces_.size(), -1);
115  neighbour_.setSize(faces_.size(), -1);
116 
117  boolList markedFaces(faces_.size(), false);
118 
119  label nInternalFaces = 0;
120 
121  forAll(c, celli)
122  {
123  // get reference to face labels for current cell
124  const labelList& cellfaces = c[celli];
125 
126  forAll(cellfaces, facei)
127  {
128  if (cellfaces[facei] < 0)
129  {
131  << "Illegal face label " << cellfaces[facei]
132  << " in cell " << celli
133  << exit(FatalError);
134  }
135 
136  if (!markedFaces[cellfaces[facei]])
137  {
138  // First visit: owner
139  owner_[cellfaces[facei]] = celli;
140  markedFaces[cellfaces[facei]] = true;
141  }
142  else
143  {
144  // Second visit: neighbour
145  neighbour_[cellfaces[facei]] = celli;
146  nInternalFaces++;
147  }
148  }
149  }
150 
151  // The neighbour array is initialised with the same length as the owner
152  // padded with -1s and here it is truncated to the correct size of the
153  // number of internal faces.
154  neighbour_.setSize(nInternalFaces);
155 
156  // Reset the primitiveMesh
158  (
159  points_.size(),
160  neighbour_.size(),
161  owner_.size(),
162  c.size(),
163  c
164  );
165 
166  string meshInfo =
167  "nPoints: " + Foam::name(nPoints())
168  + " nCells: " + Foam::name(nCells())
169  + " nFaces: " + Foam::name(nFaces())
170  + " nInternalFaces: " + Foam::name(this->nInternalFaces());
171 
172  owner_.note() = meshInfo;
173  neighbour_.note() = meshInfo;
174 }
175 
176 
177 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
string & note()
Return non-constant access to the optional note.
Definition: IOobject.H:328
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
label nCells() const
label nPoints() const
label nInternalFaces() const
label nFaces() const
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label nPoints
#define InfoInFunction
Report an information message using Foam::Info.
const dimensionedScalar c
Speed of light in a vacuum.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
List< cell > cellList
list of cells
Definition: cellList.H:42
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
error FatalError