SloanRenumber.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) 2012-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  Adapted from Boost graph/example/sloan_ordering.cpp
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "SloanRenumber.H"
30 #include "decompositionMethod.H"
31 #include "processorPolyPatch.H"
32 #include "syncTools.H"
33 
34 #include <boost/config.hpp>
35 #include <vector>
36 #include <iostream>
37 #include <boost/graph/adjacency_list.hpp>
38 #include <boost/graph/sloan_ordering.hpp>
39 #include <boost/graph/properties.hpp>
40 #include <boost/graph/bandwidth.hpp>
41 #include <boost/graph/profile.hpp>
42 #include <boost/graph/wavefront.hpp>
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 using namespace boost;
47 using namespace std;
48 
49 //Defining the graph type
50 typedef adjacency_list
51 <
52  setS,
53  vecS,
54  undirectedS,
55  property
56  <
57  vertex_color_t,
58  default_color_type,
59  property
60  <
61  vertex_degree_t,
63  property
64  <
65  vertex_priority_t,
66  Foam::scalar
67  >
68  >
69  >
71 
72 typedef graph_traits<Graph>::vertex_descriptor Vertex;
73 typedef graph_traits<Graph>::vertices_size_type size_type;
74 
75 namespace Foam
76 {
77  defineTypeNameAndDebug(SloanRenumber, 0);
78 
80  (
81  renumberMethod,
82  SloanRenumber,
83  dictionary
84  );
85 }
86 
87 
88 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
89 
90 Foam::SloanRenumber::SloanRenumber(const dictionary& renumberDict)
91 :
92  renumberMethod(renumberDict),
93  reverse_
94  (
95  renumberDict.optionalSubDict
96  (
97  typeName + "Coeffs"
98  ).lookupOrDefault<Switch>("reverse", false)
99  )
100 {}
101 
102 
103 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
104 
106 (
107  const polyMesh& mesh,
108  const pointField& points
109 ) const
110 {
111  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
112 
113  // Construct graph : faceOwner + connections across cyclics.
114 
115  // Determine neighbour cell
116  labelList nbr(mesh.nFaces()-mesh.nInternalFaces(), -1);
117  forAll(pbm, patchi)
118  {
119  if (pbm[patchi].coupled() && !isA<processorPolyPatch>(pbm[patchi]))
120  {
122  (
123  nbr,
124  pbm[patchi].size(),
125  pbm[patchi].start()-mesh.nInternalFaces()
126  ) = pbm[patchi].faceCells();
127  }
128  }
130 
131 
132  Graph G(mesh.nCells());
133 
134  // Add internal faces
135  forAll(mesh.faceNeighbour(), facei)
136  {
137  add_edge(mesh.faceOwner()[facei], mesh.faceNeighbour()[facei], G);
138  }
139  // Add cyclics
140  forAll(pbm, patchi)
141  {
142  if
143  (
144  pbm[patchi].coupled()
145  && !isA<processorPolyPatch>(pbm[patchi])
146  && refCast<const coupledPolyPatch>(pbm[patchi]).owner()
147  )
148  {
149  const labelUList& faceCells = pbm[patchi].faceCells();
150  forAll(faceCells, i)
151  {
152  label bFacei = pbm[patchi].start()+i-mesh.nInternalFaces();
153  label nbrCelli = nbr[bFacei];
154 
155  if (faceCells[i] < nbrCelli)
156  {
157  add_edge(faceCells[i], nbrCelli, G);
158  }
159  else
160  {
161  add_edge(nbrCelli, faceCells[i], G);
162  }
163  }
164  }
165  }
166 
167 
168  // Creating two iterators over the vertices
169  graph_traits<Graph>::vertex_iterator ui, ui_end;
170 
171  // Creating a property_map with the degrees of the degrees of each vertex
172  property_map<Graph,vertex_degree_t>::type deg = get(vertex_degree, G);
173  for (boost::tie(ui, ui_end) = vertices(G); ui != ui_end; ++ui)
174  deg[*ui] = degree(*ui, G);
175 
176  // Creating a property_map for the indices of a vertex
177  property_map<Graph, vertex_index_t>::type index_map = get(vertex_index, G);
178 
179  // Creating a vector of vertices
180  std::vector<Vertex> sloan_order(num_vertices(G));
181 
182  sloan_ordering
183  (
184  G,
185  sloan_order.begin(),
186  get(vertex_color, G),
187  make_degree_map(G),
188  get(vertex_priority, G)
189  );
190 
191  labelList orderedToOld(sloan_order.size());
192  forAll(orderedToOld, c)
193  {
194  orderedToOld[c] = index_map[sloan_order[c]];
195  }
196 
197  if (reverse_)
198  {
199  reverse(orderedToOld);
200  }
201 
202  return orderedToOld;
203 }
204 
205 
207 (
208  const labelListList& cellCells,
209  const pointField& points
210 ) const
211 {
212  Graph G(cellCells.size());
213 
214  forAll(cellCells, celli)
215  {
216  const labelList& nbrs = cellCells[celli];
217  forAll(nbrs, i)
218  {
219  if (nbrs[i] > celli)
220  {
221  add_edge(celli, nbrs[i], G);
222  }
223  }
224  }
225 
226  // Creating two iterators over the vertices
227  graph_traits<Graph>::vertex_iterator ui, ui_end;
228 
229  // Creating a property_map with the degrees of the degrees of each vertex
230  property_map<Graph,vertex_degree_t>::type deg = get(vertex_degree, G);
231  for (boost::tie(ui, ui_end) = vertices(G); ui != ui_end; ++ui)
232  deg[*ui] = degree(*ui, G);
233 
234  // Creating a property_map for the indices of a vertex
235  property_map<Graph, vertex_index_t>::type index_map = get(vertex_index, G);
236 
237  // Creating a vector of vertices
238  std::vector<Vertex> sloan_order(num_vertices(G));
239 
240  sloan_ordering
241  (
242  G,
243  sloan_order.begin(),
244  get(vertex_color, G),
245  make_degree_map(G),
246  get(vertex_priority, G)
247  );
248 
249  labelList orderedToOld(sloan_order.size());
250  forAll(orderedToOld, c)
251  {
252  orderedToOld[c] = index_map[sloan_order[c]];
253  }
254 
255  if (reverse_)
256  {
257  reverse(orderedToOld);
258  }
259 
260  return orderedToOld;
261 }
262 
263 
264 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
graph_traits< Graph >::vertex_descriptor Vertex
Definition: SloanRenumber.C:72
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
type
Types of root.
Definition: Roots.H:52
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
virtual labelList renumber(const pointField &) const
Return the order in which cells need to be visited, i.e.
Definition: SloanRenumber.H:87
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to hash-table of functions with typename as the key.
label nInternalFaces() const
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1047
const dimensionedScalar G
Newtonian constant of gravitation.
label nFaces() const
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
label nCells() const
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
Macros for easy insertion into run-time selection tables.
pointField vertices(const blockVertexList &bvl)
Abstract base class for renumbering.
A List obtained as a section of another List.
Definition: SubList.H:53
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:119
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1041
adjacency_list< setS, vecS, undirectedS, property< vertex_color_t, default_color_type, property< vertex_degree_t, Foam::label, property< vertex_priority_t, Foam::scalar > > >> Graph
Definition: SloanRenumber.C:70
graph_traits< Graph >::vertices_size_type size_type
Definition: SloanRenumber.C:73
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:61
Foam::polyBoundaryMesh.
void reverse(UList< T > &, const label n)
Definition: UListI.H:322
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label patchi
const dimensionedScalar c
Speed of light in a vacuum.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &l)
Swap coupled boundary face values.
Definition: syncTools.H:430
Namespace for OpenFOAM.