FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 <OpenFOAM/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  (
67  "PrimitivePatch<Face, FaceList, PointField, PointType>::"
68  "calcAddressing()"
69  ) << "addressing already calculated"
70  << abort(FatalError);
71  }
72 
73  // get reference to localFaces
74  const List<Face>& locFcs = localFaces();
75 
76  // get reference to pointFaces
77  const labelListList& pf = pointFaces();
78 
79  // Guess the max number of edges and neighbours for a face
80  label maxEdges = 0;
81  forAll (locFcs, faceI)
82  {
83  maxEdges += locFcs[faceI].size();
84  }
85 
86  // create the lists for the various results. (resized on completion)
87  edgesPtr_ = new edgeList(maxEdges);
88  edgeList& edges = *edgesPtr_;
89 
90  edgeFacesPtr_ = new labelListList(maxEdges);
91  labelListList& edgeFaces = *edgeFacesPtr_;
92 
93  // faceFaces created using a dynamic list. Cannot guess size because
94  // of multiple connections
95  List<DynamicList<label> > ff(locFcs.size());
96 
97  faceEdgesPtr_ = new labelListList(locFcs.size());
98  labelListList& faceEdges = *faceEdgesPtr_;
99 
100  // count the number of face neighbours
101  labelList noFaceFaces(locFcs.size());
102 
103  // initialise the lists of subshapes for each face to avoid duplication
104  edgeListList faceIntoEdges(locFcs.size());
105 
106  forAll (locFcs, faceI)
107  {
108  faceIntoEdges[faceI] = locFcs[faceI].edges();
109 
110  labelList& curFaceEdges = faceEdges[faceI];
111  curFaceEdges.setSize(faceIntoEdges[faceI].size());
112 
113  forAll (curFaceEdges, faceEdgeI)
114  {
115  curFaceEdges[faceEdgeI] = -1;
116  }
117  }
118 
119  // This algorithm will produce a separated list of edges, internal edges
120  // starting from 0 and boundary edges starting from the top and
121  // growing down.
122 
123  label nEdges = 0;
124 
125  bool found = false;
126 
127  // Note that faceIntoEdges is sorted acc. to local vertex numbering
128  // in face (i.e. curEdges[0] is edge between f[0] and f[1])
129 
130  // For all local faces ...
131  forAll (locFcs, faceI)
132  {
133  // Get reference to vertices of current face and corresponding edges.
134  const Face& curF = locFcs[faceI];
135  const edgeList& curEdges = faceIntoEdges[faceI];
136 
137  // Record the neighbour face. Multiple connectivity allowed
138  List<DynamicList<label> > neiFaces(curF.size());
139  List<DynamicList<label> > edgeOfNeiFace(curF.size());
140 
141  label nNeighbours = 0;
142 
143  // For all edges ...
144  forAll (curEdges, edgeI)
145  {
146  // If the edge is already detected, skip
147  if (faceEdges[faceI][edgeI] >= 0) continue;
148 
149  found = false;
150 
151  // Set reference to the current edge
152  const edge& e = curEdges[edgeI];
153 
154  // Collect neighbours for the current face vertex.
155 
156  const labelList& nbrFaces = pf[e.start()];
157 
158  forAll (nbrFaces, nbrFaceI)
159  {
160  // set reference to the current neighbour
161  label curNei = nbrFaces[nbrFaceI];
162 
163  // Reject neighbours with the lower label
164  if (curNei > faceI)
165  {
166  // get the reference to subshapes of the neighbour
167  const edgeList& searchEdges = faceIntoEdges[curNei];
168 
169  forAll (searchEdges, neiEdgeI)
170  {
171  if (searchEdges[neiEdgeI] == e)
172  {
173  // Match
174  found = true;
175 
176  neiFaces[edgeI].append(curNei);
177  edgeOfNeiFace[edgeI].append(neiEdgeI);
178 
179  // Record faceFaces both ways
180  ff[faceI].append(curNei);
181  ff[curNei].append(faceI);
182 
183  // Keep searching due to multiple connectivity
184  }
185  }
186  }
187  } // End of neighbouring faces
188 
189  if (found)
190  {
191  // Register another detected internal edge
192  nNeighbours++;
193  }
194  } // End of current edges
195 
196  // Add the edges in increasing number of neighbours.
197  // Note: for multiply connected surfaces, the lower index neighbour for
198  // an edge will come first.
199 
200  // Add the faces in the increasing order of neighbours
201  for (label neiSearch = 0; neiSearch < nNeighbours; neiSearch++)
202  {
203  // Find the lowest neighbour which is still valid
204  label nextNei = -1;
205  label minNei = locFcs.size();
206 
207  forAll (neiFaces, nfI)
208  {
209  if (neiFaces[nfI].size() && neiFaces[nfI][0] < minNei)
210  {
211  nextNei = nfI;
212  minNei = neiFaces[nfI][0];
213  }
214  }
215 
216  if (nextNei > -1)
217  {
218  // Add the face to the list of faces
219  edges[nEdges] = curEdges[nextNei];
220 
221  // Set face-edge and face-neighbour-edge to current face label
222  faceEdges[faceI][nextNei] = nEdges;
223 
224  DynamicList<label>& cnf = neiFaces[nextNei];
225  DynamicList<label>& eonf = edgeOfNeiFace[nextNei];
226 
227  // Set edge-face addressing
228  labelList& curEf = edgeFaces[nEdges];
229  curEf.setSize(cnf.size() + 1);
230  curEf[0] = faceI;
231 
232  forAll (cnf, cnfI)
233  {
234  faceEdges[cnf[cnfI]][eonf[cnfI]] = nEdges;
235 
236  curEf[cnfI + 1] = cnf[cnfI];
237  }
238 
239  // Stop the neighbour from being used again
240  cnf.clear();
241  eonf.clear();
242 
243  // Increment number of faces counter
244  nEdges++;
245  }
246  else
247  {
249  (
250  "PrimitivePatch<Face, FaceList, PointField, PointType>::"
251  "calcAddressing()"
252  ) << "Error in internal edge insertion"
253  << abort(FatalError);
254  }
255  }
256  }
257 
258  nInternalEdges_ = nEdges;
259 
260  // Do boundary faces
261 
262  forAll (faceEdges, faceI)
263  {
264  labelList& curEdges = faceEdges[faceI];
265 
266  forAll (curEdges, edgeI)
267  {
268  if (curEdges[edgeI] < 0)
269  {
270  // Grab edge and faceEdge
271  edges[nEdges] = faceIntoEdges[faceI][edgeI];
272  curEdges[edgeI] = nEdges;
273 
274  // Add edgeFace
275  labelList& curEf = edgeFaces[nEdges];
276  curEf.setSize(1);
277  curEf[0] = faceI;
278 
279  nEdges++;
280  }
281  }
282  }
283 
284  // edges
285  edges.setSize(nEdges);
286 
287  // edgeFaces list
288  edgeFaces.setSize(nEdges);
289 
290  // faceFaces list
291  faceFacesPtr_ = new labelListList(locFcs.size());
292  labelListList& faceFaces = *faceFacesPtr_;
293 
294  forAll (faceFaces, faceI)
295  {
296  faceFaces[faceI].transfer(ff[faceI]);
297  }
298 
299 
300  if (debug)
301  {
302  Info<< "PrimitivePatch<Face, FaceList, PointField, PointType>::"
303  << "calcAddressing() : finished calculating patch addressing"
304  << endl;
305  }
306 }
307 
308 
309 // ************************ vim: set sw=4 sts=4 et: ************************ //