structuredRenumber.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) 2018-2020 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\*---------------------------------------------------------------------------*/
28
29#include "structuredRenumber.H"
31#include "topoDistanceData.H"
32#include "fvMeshSubset.H"
34
35// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36
37namespace Foam
38{
40
42 (
46 );
47}
48
49
50// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51
53(
54 const dictionary& dict
55)
56:
58 coeffsDict_(dict.optionalSubDict(typeName + "Coeffs")),
59 patches_(coeffsDict_.get<wordRes>("patches")),
60 nLayers_(coeffsDict_.getOrDefault<label>("nLayers", labelMax)),
61 depthFirst_(coeffsDict_.get<bool>("depthFirst")),
62 reverse_(coeffsDict_.get<bool>("reverse")),
63 method_(renumberMethod::New(coeffsDict_))
64{}
65
66
67// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
68
70(
71 const label a,
72 const label b
73)
74{
75 const topoDistanceData<label>& ta = distance_[a];
76 const topoDistanceData<label>& tb = distance_[b];
77
78 int dummy;
79
80 if (ta.valid(dummy))
81 {
82 if (tb.valid(dummy))
83 {
84 if (depthFirst_)
85 {
86 if (ta.data() < tb.data())
87 {
88 // Sort column first
89 return true;
90 }
91 else if (ta.data() > tb.data())
92 {
93 return false;
94 }
95 else
96 {
97 // Same column. Sort according to layer
98 return ta.distance() < tb.distance();
99 }
100 }
101 else
102 {
103 if (ta.distance() < tb.distance())
104 {
105 return true;
106 }
107 else if (ta.distance() > tb.distance())
108 {
109 return false;
110 }
111 else
112 {
113 // Same layer; sort according to current values
114 return ta.data() < tb.data();
115 }
116 }
117 }
118 else
119 {
120 return true;
121 }
122 }
123 else
124 {
125 if (tb.valid(dummy))
126 {
127 return false;
128 }
129 else
130 {
131 // Both not valid; fall back to cell indices for sorting
132 return order_[a] < order_[b];
133 }
134 }
135}
136
137
139(
140 const polyMesh& mesh,
141 const pointField& points
142) const
143{
144 if (points.size() != mesh.nCells())
145 {
147 << "Number of points " << points.size()
148 << " should equal the number of cells " << mesh.nCells()
149 << exit(FatalError);
150 }
151
152 const polyBoundaryMesh& pbm = mesh.boundaryMesh();
153 const labelHashSet patchIDs(pbm.patchSet(patches_));
154
155 label nFaces = 0;
156 for (const label patchi : patchIDs)
157 {
158 nFaces += pbm[patchi].size();
159 }
160
161
162 // Extract a submesh.
163 labelHashSet patchCells(2*nFaces);
164 for (const label patchId : patchIDs)
165 {
166 patchCells.insert(pbm[patchId].faceCells());
167 }
168
169 label nTotalSeeds = returnReduce(patchCells.size(), sumOp<label>());
170
171 label nTotalCells = mesh.globalData().nTotalCells();
172 const label nLayers = nTotalCells/nTotalSeeds;
173
174 Info<< type() << " : seeding " << nTotalSeeds
175 << " cells on (estimated) " << nLayers << " layers" << nl
176 << endl;
177
178
179 // Work array. Used here to temporarily store the original-to-ordered
180 // index. Later on used to store the ordered-to-original.
181 labelList orderedToOld(mesh.nCells(), -1);
182
183 // Subset the layer of cells next to the patch
184 {
185 fvMeshSubset subsetter
186 (
187 dynamic_cast<const fvMesh&>(mesh),
188 patchCells
189 );
190 const fvMesh& subMesh = subsetter.subMesh();
191
192 pointField subPoints(points, subsetter.cellMap());
193
194 // Locally renumber the layer of cells
195 labelList subOrder(method_().renumber(subMesh, subPoints));
196
197 labelList subOrigToOrdered(invert(subOrder.size(), subOrder));
198
199 globalIndex globalSubCells(subOrder.size());
200
201 // Transfer to final decomposition and convert into global numbering
202 forAll(subOrder, i)
203 {
204 orderedToOld[subsetter.cellMap()[i]] =
205 globalSubCells.toGlobal(subOrigToOrdered[i]);
206 }
207 }
208
209
210 // Walk sub-ordering (=column index) out.
211 labelList patchFaces(nFaces);
213 nFaces = 0;
214 for (const label patchi : patchIDs)
215 {
216 const polyPatch& pp = pbm[patchi];
217 const labelUList& fc = pp.faceCells();
218 forAll(fc, i)
219 {
220 patchFaces[nFaces] = pp.start()+i;
222 (
223 0, // distance: layer
224 orderedToOld[fc[i]] // passive data: global column
225 );
226 nFaces++;
227 }
228 }
229
230 // Field on cells and faces.
233
234 // Propagate information inwards
236 (
237 mesh,
238 patchFaces,
239 patchData,
240 faceData,
241 cellData,
242 0
243 );
244
245 deltaCalc.iterate(nLayers_);
246
247 Info<< type() << " : did not visit "
248 << deltaCalc.nUnvisitedCells()
249 << " cells out of " << nTotalCells
250 << "; using " << method_().type() << " renumbering for these" << endl;
251
252 // Get cell order using the method(). These values will get overwritten
253 // by any visited cell so are used only if the number of nLayers is limited.
254 labelList oldToOrdered
255 (
256 invert
257 (
258 mesh.nCells(),
259 method_().renumber(mesh, points)
260 )
261 );
262
263 // Use specialised sorting to sorted either layers or columns first
264 // Done so that at no point we need to combine both into a single
265 // index and we might run out of label size.
267 (
268 cellData,
269 orderedToOld,
270 layerLess(depthFirst_, oldToOrdered, cellData)
271 );
272
273 // Return furthest away cell first
274 if (reverse_)
275 {
276 reverse(orderedToOld);
277 }
278
279 return orderedToOld;
280}
281
282
283// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
label nUnvisitedCells() const
virtual label iterate(const label maxIter)
Iterate until no changes or maxIter reached.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition: HashSet.H:191
label size() const noexcept
The number of elements in table.
Definition: HashTableI.H:52
Version of FaceCellWave that walks through prismatic cells only.
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:106
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
friend Ostream & operator(Ostream &, const faMatrix< Type > &)
Smooth ATC in cells next to a set of patches supplied by type.
Definition: faceCells.H:59
Holds a reference to the original mesh (the baseMesh) and optionally to a subset of that mesh (the su...
Definition: fvMeshSubset.H:80
const labelList & cellMap() const
Return cell map.
Definition: fvMeshSubsetI.H:91
const fvMesh & subMesh() const
Return reference to subset mesh.
Definition: fvMeshSubsetI.H:48
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:68
label toGlobal(const label i) const
From local to global index.
Definition: globalIndexI.H:260
label nTotalCells() const noexcept
Return total number of cells in decomposed mesh.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool useGroups=true) const
Return the set of patch IDs corresponding to the given names.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1310
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:456
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:364
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:371
label nCells() const noexcept
Number of mesh cells.
label nFaces() const noexcept
Number of mesh faces.
Abstract base class for renumbering.
Less function class that can be used for sorting according to.
Renumbering according to mesh layers. depthFirst = true: first column gets ids 0.....
virtual labelList renumber(const pointField &) const
For use with FaceCellWave. Determines topological distance to starting faces. Templated on passive tr...
bool valid(TrackingData &td) const
Changed or contains original (invalid) value.
const Type & data() const
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:54
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
bool
Definition: EEqn.H:20
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
const labelList nFaces(UPstream::listGatherValues< label >(aMesh.nFaces()))
const pointField & points
label patchId(-1)
Namespace for OpenFOAM.
void reverse(UList< T > &list, const label n)
Reverse the first n elements of the list.
Definition: UListI.H:449
messageStream Info
Information stream (stdout output on master, null elsewhere)
IntListType renumber(const labelUList &oldToNew, const IntListType &input)
Renumber the values (not the indices) of a list.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
constexpr label labelMax
Definition: label.H:61
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
error FatalError
labelList invert(const label len, const labelUList &map)
Create an inverse one-to-one mapping.
Definition: ListOps.C:36
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
volScalarField & b
Definition: createFields.H:27
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333