cellShapeI.H
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 "Istream.H"
27 #include "cell.H"
28 #include "cellModeller.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
33 :
34  m(nullptr)
35 {}
36 
37 
39 (
40  const cellModel& M,
41  const labelList& l,
42  const bool doCollapse
43 )
44 :
45  labelList(l),
46  m(&M)
47 {
48  if (doCollapse)
49  {
50  collapse();
51  }
52 }
53 
54 
56 (
57  const word& model,
58  const labelList& l,
59  const bool doCollapse
60 )
61 :
62  labelList(l),
63  m(cellModeller::lookup(model))
64 {
65  if (doCollapse)
66  {
67  collapse();
68  }
69 }
70 
71 
73 {
74  is >> *this;
75 }
76 
77 
79 {
80  return autoPtr<cellShape>(new cellShape(*this));
81 }
82 
83 
84 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
85 
87 (
88  const pointField& meshPoints
89 ) const
90 {
91  // There are as many points as there labels for them
92  pointField p(size());
93 
94  // For each point in list, set it to the point in 'pnts' addressed
95  // by 'labs'
96  forAll(p, i)
97  {
98  p[i] = meshPoints[operator[](i)];
99  }
100 
101  // Return list
102  return p;
103 }
104 
105 
107 {
108  return *m;
109 }
110 
111 
113 (
114  const faceList& allFaces,
115  const cell& cFaces
116 ) const
117 {
118  // Faces in model order
119  faceList localFaces(faces());
120 
121  // Do linear match (usually cell shape is low complexity)
122 
123  labelList modelToMesh(localFaces.size(), -1);
124 
125  forAll(localFaces, i)
126  {
127  const face& localF = localFaces[i];
128 
129  forAll(cFaces, j)
130  {
131  label meshFacei = cFaces[j];
132 
133  if (allFaces[meshFacei] == localF)
134  {
135  modelToMesh[i] = meshFacei;
136 
137  break;
138  }
139  }
140  }
141 
142  return modelToMesh;
143 }
144 
145 
147 (
148  const edgeList& allEdges,
149  const labelList& cEdges
150 ) const
151 {
152  // Edges in model order
153  edgeList localEdges(edges());
154 
155  // Do linear match (usually cell shape is low complexity)
156 
157  labelList modelToMesh(localEdges.size(), -1);
158 
159  forAll(localEdges, i)
160  {
161  const edge& e = localEdges[i];
162 
163  forAll(cEdges, j)
164  {
165  label edgeI = cEdges[j];
166 
167  if (allEdges[edgeI] == e)
168  {
169  modelToMesh[i] = edgeI;
170 
171  break;
172  }
173  }
174  }
175 
176  return modelToMesh;
177 }
178 
179 
181 {
182  return m->faces(*this);
183 }
184 
185 
187 {
188  faceList oldFaces(faces());
189 
190  faceList newFaces(oldFaces.size());
191  label newFacei = 0;
192 
193  forAll(oldFaces, oldFacei)
194  {
195  const face& f = oldFaces[oldFacei];
196 
197  face& newF = newFaces[newFacei];
198 
199  newF.setSize(f.size());
200 
201  label newFp = 0;
202  label prevVertI = -1;
203 
204  forAll(f, fp)
205  {
206  label vertI = f[fp];
207 
208  if (vertI != prevVertI)
209  {
210  newF[newFp++] = vertI;
211 
212  prevVertI = vertI;
213  }
214  }
215 
216  if ((newFp > 1) && (newF[newFp-1] == newF[0]))
217  {
218  --newFp;
219  }
220 
221  if (newFp > 2)
222  {
223  // Size face and go to next one
224  newF.setSize(newFp);
225 
226  newFacei++;
227  }
228  }
229  newFaces.setSize(newFacei);
230 
231  return newFaces;
232 }
233 
234 
236 {
237  return m->nFaces();
238 }
239 
240 
242 {
243  return m->edges(*this);
244 }
245 
246 
248 {
249  return m->nEdges();
250 }
251 
252 
254 {
255  return size();
256 }
257 
258 
260 {
261  return m->centre(*this, points);
262 }
263 
264 
265 inline Foam::scalar Foam::cellShape::mag(const pointField& points) const
266 {
267  return m->mag(*this, points);
268 }
269 
270 
271 // ************************************************************************* //
scalar mag(const pointField &) const
Scalar magnitude.
Definition: cellShapeI.H:265
label nFaces() const
Return number of faces.
Definition: cellModelI.H:62
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
static const cellModel * lookup(const word &)
Look up a model by name and return a pointer to the model or nullptr.
Definition: cellModeller.C:100
labelList meshFaces(const faceList &allFaces, const cell &) const
Mesh face labels of this cell (in order of model)
Definition: cellShapeI.H:113
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
T & operator[](const label)
Return element of UList.
Definition: UListI.H:167
void collapse()
Collapse shape to correct one after removing duplicate vertices.
Definition: cellShape.C:31
label nEdges() const
Return number of edges.
Definition: cellModelI.H:56
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
vector centre(const labelList &pointLabels, const pointField &points) const
Vector centroid.
Definition: cellModel.C:32
autoPtr< cellShape > clone() const
Clone.
Definition: cellShapeI.H:78
pointField points(const pointField &meshPoints) const
Return the points corresponding to this cellShape.
Definition: cellShapeI.H:87
labelList meshEdges(const edgeList &allEdges, const labelList &) const
Mesh edge labels of this cell (in order of model)
Definition: cellShapeI.H:147
faceList faces() const
Faces of this cell.
Definition: cellShapeI.H:180
cellShape()
Construct null.
Definition: cellShapeI.H:32
label nEdges() const
Number of edges.
Definition: cellShapeI.H:247
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
A class for handling words, derived from string.
Definition: word.H:59
List< label > labelList
A List of labels.
Definition: labelList.H:56
faceList collapsedFaces() const
Collapsed faces of this cell.
Definition: cellShapeI.H:186
labelList f(nPoints)
edgeList edges() const
Edges of this cellShape.
Definition: cellShapeI.H:241
void setSize(const label)
Reset size of List.
Definition: List.C:281
const cellModel & model() const
Model reference.
Definition: cellShapeI.H:106
scalar mag(const labelList &pointLabels, const pointField &points) const
Cell volume.
Definition: cellModel.C:89
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:57
point centre(const pointField &) const
Centroid of the cell.
Definition: cellShapeI.H:259
Maps a geometry to a set of cell primitives, which enables geometric cell data to be calculated witho...
Definition: cellModel.H:64
edgeList edges(const labelList &pointLabels) const
Return list of edges.
Definition: cellModelI.H:70
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
label nFaces() const
Number of faces.
Definition: cellShapeI.H:235
volScalarField & p
label size() const
Return the number of elements in the UList.
Definition: ListI.H:171
label nPoints() const
Number of points.
Definition: cellShapeI.H:253
faceList faces(const labelList &pointLabels) const
Return list of faces.
Definition: cellModelI.H:97