37#include <boost/config.hpp>
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>
74typedef graph_traits<Graph>::vertex_descriptor
Vertex;
75typedef graph_traits<Graph>::vertices_size_type
size_type;
97 dict.optionalSubDict(typeName +
"Coeffs")
98 .getOrDefault(
"reverse", false)
110 using namespace Foam;
113 graph_traits<Graph>::vertex_iterator ui, ui_end;
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)
119 deg[*ui] = degree(*ui, G);
123 property_map<Graph, vertex_index_t>::type index_map = get(vertex_index, G);
126 std::vector<Vertex> sloan_order(num_vertices(G));
132 get(vertex_color, G),
134 get(vertex_priority, G)
137 labelList orderedToOld(sloan_order.size());
140 orderedToOld[c] = index_map[sloan_order[c]];
168 if (pp.coupled() && !isA<processorPolyPatch>(pp))
190 && !isA<processorPolyPatch>(pp)
191 && refCast<const coupledPolyPatch>(pp).owner()
194 label bFacei = pp.offset();
196 for (
const label ownCelli : pp.faceCells())
198 const label nbrCelli = nbr[bFacei];
201 if (ownCelli < nbrCelli)
203 add_edge(ownCelli, nbrCelli, G);
207 add_edge(nbrCelli, ownCelli, G);
213 return renumberImpl(G, reverse_);
227 const auto& neighbours = cellCells[celli];
229 for (
const label nbr : neighbours)
233 add_edge(celli, nbr, G);
238 return renumberImpl(G, reverse_);
252 const auto& neighbours = cellCells[celli];
254 for (
const label nbr : neighbours)
258 add_edge(celli, nbr, G);
263 return renumberImpl(G, reverse_);
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
graph_traits< Graph >::vertices_size_type size_type
graph_traits< Graph >::vertex_descriptor Vertex
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.
virtual labelList renumber(const pointField &) const
A List obtained as a section of another List.
void size(const label n)
Older name for setAddressableSize.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh consisting of general polyhedral cells.
virtual const labelList & faceOwner() const
Return face owner.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const labelList & faceNeighbour() const
Return face neighbour.
A patch is a list of labels that address the faces in the global face list.
label nBoundaryFaces() const noexcept
Number of boundary faces (== nFaces - nInternalFaces)
label nCells() const noexcept
Number of mesh cells.
Abstract base class for renumbering.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
void reverse(UList< T > &list, const label n)
Reverse the first n elements of the list.
pointField vertices(const blockVertexList &bvl)
#define forAll(list, i)
Loop across all elements in list.