PrimitivePatchAddressing.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 Description
25  This function calculates the list of patch edges, defined on the list of
26  points supporting the patch. The edges are ordered:
27  - 0..nInternalEdges-1 : upper triangular order
28  - nInternalEdges.. : boundary edges (no particular order)
29 
30  Other patch addressing information is also calculated:
31  - faceFaces with neighbour faces in ascending order
32  - edgeFaces with ascending face order
33  - faceEdges sorted according to edges of a face
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #include "PrimitivePatch.H"
38 #include "DynamicList.H"
39 
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
43 template
44 <
45  class Face,
46  template<class> class FaceList,
47  class PointField,
48  class PointType
49 >
50 void
52 calcAddressing() const
53 {
54  if (debug)
55  {
56  Info<< "PrimitivePatch<Face, FaceList, PointField, PointType>::"
57  << "calcAddressing() : calculating patch addressing"
58  << endl;
59  }
60 
61  if (edgesPtr_ || faceFacesPtr_ || edgeFacesPtr_ || faceEdgesPtr_)
62  {
63  // it is considered an error to attempt to recalculate
64  // if already allocated
66  << "addressing already calculated"
67  << abort(FatalError);
68  }
69 
70  // get reference to localFaces
71  const List<Face>& locFcs = localFaces();
72 
73  // get reference to pointFaces
74  const labelListList& pf = pointFaces();
75 
76  // Guess the max number of edges and neighbours for a face
77  label maxEdges = 0;
78  forAll(locFcs, facei)
79  {
80  maxEdges += locFcs[facei].size();
81  }
82 
83  // create the lists for the various results. (resized on completion)
84  edgesPtr_ = new edgeList(maxEdges);
85  edgeList& edges = *edgesPtr_;
86 
87  edgeFacesPtr_ = new labelListList(maxEdges);
88  labelListList& edgeFaces = *edgeFacesPtr_;
89 
90  // faceFaces created using a dynamic list. Cannot guess size because
91  // of multiple connections
92  List<DynamicList<label>> ff(locFcs.size());
93 
94  faceEdgesPtr_ = new labelListList(locFcs.size());
95  labelListList& faceEdges = *faceEdgesPtr_;
96 
97  // count the number of face neighbours
98  labelList noFaceFaces(locFcs.size());
99 
100  // initialise the lists of subshapes for each face to avoid duplication
101  edgeListList faceIntoEdges(locFcs.size());
102 
103  forAll(locFcs, facei)
104  {
105  faceIntoEdges[facei] = locFcs[facei].edges();
106 
107  labelList& curFaceEdges = faceEdges[facei];
108  curFaceEdges.setSize(faceIntoEdges[facei].size());
109 
110  forAll(curFaceEdges, faceEdgeI)
111  {
112  curFaceEdges[faceEdgeI] = -1;
113  }
114  }
115 
116  // This algorithm will produce a separated list of edges, internal edges
117  // starting from 0 and boundary edges starting from the top and
118  // growing down.
119 
120  label nEdges = 0;
121 
122  bool found = false;
123 
124  // Note that faceIntoEdges is sorted acc. to local vertex numbering
125  // in face (i.e. curEdges[0] is edge between f[0] and f[1])
126 
127  // For all local faces ...
128  forAll(locFcs, facei)
129  {
130  // Get reference to vertices of current face and corresponding edges.
131  const Face& curF = locFcs[facei];
132  const edgeList& curEdges = faceIntoEdges[facei];
133 
134  // Record the neighbour face. Multiple connectivity allowed
135  List<DynamicList<label>> neiFaces(curF.size());
136  List<DynamicList<label>> edgeOfNeiFace(curF.size());
137 
138  label nNeighbours = 0;
139 
140  // For all edges ...
141  forAll(curEdges, edgeI)
142  {
143  // If the edge is already detected, skip
144  if (faceEdges[facei][edgeI] >= 0) continue;
145 
146  found = false;
147 
148  // Set reference to the current edge
149  const edge& e = curEdges[edgeI];
150 
151  // Collect neighbours for the current face vertex.
152 
153  const labelList& nbrFaces = pf[e.start()];
154 
155  forAll(nbrFaces, nbrFacei)
156  {
157  // set reference to the current neighbour
158  label curNei = nbrFaces[nbrFacei];
159 
160  // Reject neighbours with the lower label
161  if (curNei > facei)
162  {
163  // get the reference to subshapes of the neighbour
164  const edgeList& searchEdges = faceIntoEdges[curNei];
165 
166  forAll(searchEdges, neiEdgeI)
167  {
168  if (searchEdges[neiEdgeI] == e)
169  {
170  // Match
171  found = true;
172 
173  neiFaces[edgeI].append(curNei);
174  edgeOfNeiFace[edgeI].append(neiEdgeI);
175 
176  // Record faceFaces both ways
177  ff[facei].append(curNei);
178  ff[curNei].append(facei);
179 
180  // Keep searching due to multiple connectivity
181  }
182  }
183  }
184  } // End of neighbouring faces
185 
186  if (found)
187  {
188  // Register another detected internal edge
189  nNeighbours++;
190  }
191  } // End of current edges
192 
193  // Add the edges in increasing number of neighbours.
194  // Note: for multiply connected surfaces, the lower index neighbour for
195  // an edge will come first.
196 
197  // Add the faces in the increasing order of neighbours
198  for (label neiSearch = 0; neiSearch < nNeighbours; neiSearch++)
199  {
200  // Find the lowest neighbour which is still valid
201  label nextNei = -1;
202  label minNei = locFcs.size();
203 
204  forAll(neiFaces, nfI)
205  {
206  if (neiFaces[nfI].size() && neiFaces[nfI][0] < minNei)
207  {
208  nextNei = nfI;
209  minNei = neiFaces[nfI][0];
210  }
211  }
212 
213  if (nextNei > -1)
214  {
215  // Add the face to the list of faces
216  edges[nEdges] = curEdges[nextNei];
217 
218  // Set face-edge and face-neighbour-edge to current face label
219  faceEdges[facei][nextNei] = nEdges;
220 
221  DynamicList<label>& cnf = neiFaces[nextNei];
222  DynamicList<label>& eonf = edgeOfNeiFace[nextNei];
223 
224  // Set edge-face addressing
225  labelList& curEf = edgeFaces[nEdges];
226  curEf.setSize(cnf.size() + 1);
227  curEf[0] = facei;
228 
229  forAll(cnf, cnfI)
230  {
231  faceEdges[cnf[cnfI]][eonf[cnfI]] = nEdges;
232 
233  curEf[cnfI + 1] = cnf[cnfI];
234  }
235 
236  // Stop the neighbour from being used again
237  cnf.clear();
238  eonf.clear();
239 
240  // Increment number of faces counter
241  nEdges++;
242  }
243  else
244  {
246  << "Error in internal edge insertion"
247  << abort(FatalError);
248  }
249  }
250  }
251 
252  nInternalEdges_ = nEdges;
253 
254  // Do boundary faces
255 
256  forAll(faceEdges, facei)
257  {
258  labelList& curEdges = faceEdges[facei];
259 
260  forAll(curEdges, edgeI)
261  {
262  if (curEdges[edgeI] < 0)
263  {
264  // Grab edge and faceEdge
265  edges[nEdges] = faceIntoEdges[facei][edgeI];
266  curEdges[edgeI] = nEdges;
267 
268  // Add edgeFace
269  labelList& curEf = edgeFaces[nEdges];
270  curEf.setSize(1);
271  curEf[0] = facei;
272 
273  nEdges++;
274  }
275  }
276  }
277 
278  // edges
279  edges.setSize(nEdges);
280 
281  // edgeFaces list
282  edgeFaces.setSize(nEdges);
283 
284  // faceFaces list
285  faceFacesPtr_ = new labelListList(locFcs.size());
286  labelListList& faceFaces = *faceFacesPtr_;
287 
288  forAll(faceFaces, facei)
289  {
290  faceFaces[facei].transfer(ff[facei]);
291  }
292 
293 
294  if (debug)
295  {
296  Info<< "PrimitivePatch<Face, FaceList, PointField, PointType>::"
297  << "calcAddressing() : finished calculating patch addressing"
298  << endl;
299  }
300 }
301 
302 
303 // ************************************************************************* //
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
#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
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
List< edgeList > edgeListList
Definition: edgeList.H:39
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A list of faces which address into the list of points.
List< edge > edgeList
Definition: edgeList.H:38
List< label > labelList
A List of labels.
Definition: labelList.H:56
errorManip< error > abort(error &err)
Definition: errorManip.H:131
void setSize(const label)
Reset size of List.
Definition: List.C:281
messageStream Info
const FieldField< fvPatchField, Type > & ff(const FieldField< fvPatchField, Type > &bf)