SloanRenumber.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2012-2017 OpenFOAM Foundation
9 Copyright (C) 2020-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27 Adapted from Boost graph/example/sloan_ordering.cpp
28
29\*---------------------------------------------------------------------------*/
30
31#include "SloanRenumber.H"
33#include "decompositionMethod.H"
34#include "processorPolyPatch.H"
35#include "syncTools.H"
36
37#include <boost/config.hpp>
38#include <vector>
39#include <iostream>
40#include <boost/graph/adjacency_list.hpp>
41#include <boost/graph/sloan_ordering.hpp>
42#include <boost/graph/properties.hpp>
43#include <boost/graph/bandwidth.hpp>
44#include <boost/graph/profile.hpp>
45#include <boost/graph/wavefront.hpp>
46
47// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48
49using namespace boost;
50
51//Defining the graph type
52typedef adjacency_list
53<
54 setS,
55 vecS,
56 undirectedS,
57 property
58 <
59 vertex_color_t,
60 default_color_type,
61 property
62 <
63 vertex_degree_t,
64 Foam::label,
65 property
66 <
67 vertex_priority_t,
68 Foam::scalar
69 >
70 >
71 >
73
74typedef graph_traits<Graph>::vertex_descriptor Vertex;
75typedef graph_traits<Graph>::vertices_size_type size_type;
76
77namespace Foam
78{
80
82 (
86 );
87}
88
89
90// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
91
93:
95 reverse_
96 (
97 dict.optionalSubDict(typeName + "Coeffs")
98 .getOrDefault("reverse", false)
99 )
100{}
101
102
103// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
104
105namespace
106{
107
108Foam::labelList renumberImpl(Graph& G, const bool useReverse)
109{
110 using namespace Foam;
111
112 //Creating two iterators over the vertices
113 graph_traits<Graph>::vertex_iterator ui, ui_end;
114
115 //Creating a property_map with the degrees of the degrees of each vertex
116 property_map<Graph,vertex_degree_t>::type deg = get(vertex_degree, G);
117 for (boost::tie(ui, ui_end) = vertices(G); ui != ui_end; ++ui)
118 {
119 deg[*ui] = degree(*ui, G);
120 }
121
122 //Creating a property_map for the indices of a vertex
123 property_map<Graph, vertex_index_t>::type index_map = get(vertex_index, G);
124
125 //Creating a vector of vertices
126 std::vector<Vertex> sloan_order(num_vertices(G));
127
128 sloan_ordering
129 (
130 G,
131 sloan_order.begin(),
132 get(vertex_color, G),
133 make_degree_map(G),
134 get(vertex_priority, G)
135 );
136
137 labelList orderedToOld(sloan_order.size());
138 forAll(orderedToOld, c)
139 {
140 orderedToOld[c] = index_map[sloan_order[c]];
141 }
142
143 if (useReverse)
144 {
145 Foam::reverse(orderedToOld);
146 }
147
148 return orderedToOld;
149}
150
151} // End anonymous namespace
152
153
154// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
155
157(
158 const polyMesh& mesh,
159 const pointField& points
160) const
161{
162 // Construct graph : faceOwner + connections across cyclics.
163
164 // Determine neighbour cell
165 labelList nbr(mesh.nBoundaryFaces(), -1);
166 for (const polyPatch& pp : mesh.boundaryMesh())
167 {
168 if (pp.coupled() && !isA<processorPolyPatch>(pp))
169 {
170 SubList<label>(nbr, pp.size(), pp.offset()) = pp.faceCells();
171 }
172 }
174
175
176 Graph G(mesh.nCells());
177
178 // Add internal faces
179 forAll(mesh.faceNeighbour(), facei)
180 {
181 add_edge(mesh.faceOwner()[facei], mesh.faceNeighbour()[facei], G);
182 }
183
184 // Add cyclics
185 for (const polyPatch& pp : mesh.boundaryMesh())
186 {
187 if
188 (
189 pp.coupled()
190 && !isA<processorPolyPatch>(pp)
191 && refCast<const coupledPolyPatch>(pp).owner()
192 )
193 {
194 label bFacei = pp.offset();
195
196 for (const label ownCelli : pp.faceCells())
197 {
198 const label nbrCelli = nbr[bFacei];
199 ++bFacei;
200
201 if (ownCelli < nbrCelli)
202 {
203 add_edge(ownCelli, nbrCelli, G);
204 }
205 else
206 {
207 add_edge(nbrCelli, ownCelli, G);
208 }
209 }
210 }
211 }
212
213 return renumberImpl(G, reverse_);
214}
215
216
218(
219 const CompactListList<label>& cellCells,
220 const pointField&
221) const
222{
223 Graph G(cellCells.size());
224
225 forAll(cellCells, celli)
226 {
227 const auto& neighbours = cellCells[celli];
228
229 for (const label nbr : neighbours)
230 {
231 if (celli < nbr)
232 {
233 add_edge(celli, nbr, G);
234 }
235 }
236 }
237
238 return renumberImpl(G, reverse_);
239}
240
241
243(
244 const labelListList& cellCells,
245 const pointField&
246) const
247{
248 Graph G(cellCells.size());
249
250 forAll(cellCells, celli)
251 {
252 const auto& neighbours = cellCells[celli];
253
254 for (const label nbr : neighbours)
255 {
256 if (celli < nbr)
257 {
258 add_edge(celli, nbr, G);
259 }
260 }
261 }
262
263 return renumberImpl(G, reverse_);
264}
265
266
267// ************************************************************************* //
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:72
graph_traits< Graph >::vertices_size_type size_type
Definition: SloanRenumber.C:75
graph_traits< Graph >::vertex_descriptor Vertex
Definition: SloanRenumber.C:74
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A packed storage unstructured matrix of objects of type <T> using an offset table for access.
label size() const noexcept
The primary size (the number of rows/sublists)
Sloan renumbering algorithm.
Definition: SloanRenumber.H:54
virtual labelList renumber(const pointField &) const
Definition: SloanRenumber.H:90
A List obtained as a section of another List.
Definition: SubList.H:70
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1121
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:456
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1127
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
label nBoundaryFaces() const noexcept
Number of boundary faces (== nFaces - nInternalFaces)
label nCells() const noexcept
Number of mesh cells.
Abstract base class for renumbering.
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &faceValues)
Swap coupled boundary face values. Uses eqOp.
Definition: syncTools.H:445
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
const pointField & points
Namespace for OpenFOAM.
void reverse(UList< T > &list, const label n)
Reverse the first n elements of the list.
Definition: UListI.H:449
pointField vertices(const blockVertexList &bvl)
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333