structuredRenumber.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-2022 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 \*---------------------------------------------------------------------------*/
25 
26 #include "structuredRenumber.H"
28 #include "topoDistanceData.H"
29 #include "fvMeshSubset.H"
30 #include "OppositeFaceCellWave.H"
31 #include "globalIndex.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
38 
40  (
44  );
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
51 (
52  const dictionary& renumberDict
53 )
54 :
55  renumberMethod(renumberDict),
56  methodDict_(renumberDict.optionalSubDict(typeName + "Coeffs")),
57  patches_(methodDict_.lookup("patches")),
58  nLayers_(methodDict_.lookupOrDefault<label>("nLayers", labelMax)),
59  depthFirst_(methodDict_.lookup("depthFirst")),
60  method_(renumberMethod::New(methodDict_)),
61  reverse_(methodDict_.lookup("reverse"))
62 {}
63 
64 
65 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
66 
67 bool Foam::structuredRenumber::layerLess::operator()
68 (
69  const label a,
70  const label b
71 )
72 {
73  const topoDistanceData& ta = distance_[a];
74  const topoDistanceData& tb = distance_[b];
75 
76  int dummy;
77 
78  if (ta.valid(dummy))
79  {
80  if (tb.valid(dummy))
81  {
82  if (depthFirst_)
83  {
84  if (ta.data() < tb.data())
85  {
86  // Sort column first
87  return true;
88  }
89  else if (ta.data() > tb.data())
90  {
91  return false;
92  }
93  else
94  {
95  // Same column. Sort according to layer
96  return ta.distance() < tb.distance();
97  }
98  }
99  else
100  {
101  if (ta.distance() < tb.distance())
102  {
103  return true;
104  }
105  else if (ta.distance() > tb.distance())
106  {
107  return false;
108  }
109  else
110  {
111  // Same layer; sort according to current values
112  return ta.data() < tb.data();
113  }
114  }
115  }
116  else
117  {
118  return true;
119  }
120  }
121  else
122  {
123  if (tb.valid(dummy))
124  {
125  return false;
126  }
127  else
128  {
129  // Both not valid; fall back to cell indices for sorting
130  return order_[a] < order_[b];
131  }
132  }
133 }
134 
135 
137 (
138  const polyMesh& mesh,
139  const pointField& points
140 ) const
141 {
142  if (points.size() != mesh.nCells())
143  {
145  << "Number of points " << points.size()
146  << " should equal the number of cells " << mesh.nCells()
147  << exit(FatalError);
148  }
149 
150  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
151  const labelHashSet patchIDs(pbm.patchSet(patches_));
152 
153  label nFaces = 0;
154  forAllConstIter(labelHashSet, patchIDs, iter)
155  {
156  nFaces += pbm[iter.key()].size();
157  }
158 
159 
160  // Extract a submesh.
161  labelHashSet patchCells(2*nFaces);
162  forAllConstIter(labelHashSet, patchIDs, iter)
163  {
164  const labelUList& fc = pbm[iter.key()].faceCells();
165  forAll(fc, i)
166  {
167  patchCells.insert(fc[i]);
168  }
169  }
170 
171  label nTotalSeeds = returnReduce(patchCells.size(), sumOp<label>());
172 
173  label nTotalCells = mesh.globalData().nTotalCells();
174  const label nLayers = nTotalCells/nTotalSeeds;
175 
176  Info<< type() << " : seeding " << nTotalSeeds
177  << " cells on (estimated) " << nLayers << " layers" << nl
178  << endl;
179 
180 
181  // Work array. Used here to temporarily store the original-to-ordered
182  // index. Later on used to store the ordered-to-original.
183  labelList orderedToOld(mesh.nCells(), -1);
184 
185  // Subset the layer of cells next to the patch
186  {
187  fvMeshSubset subsetter(dynamic_cast<const fvMesh&>(mesh));
188  subsetter.setLargeCellSubset(patchCells);
189  const fvMesh& subMesh = subsetter.subMesh();
190 
191  pointField subPoints(points, subsetter.cellMap());
192 
193  // Locally renumber the layer of cells
194  labelList subOrder(method_().renumber(subMesh, subPoints));
195 
196  labelList subOrigToOrdered(invert(subOrder.size(), subOrder));
197 
198  globalIndex globalSubCells(subOrder.size());
199 
200  // Transfer to final decomposition and convert into global numbering
201  forAll(subOrder, i)
202  {
203  orderedToOld[subsetter.cellMap()[i]] =
204  globalSubCells.toGlobal(subOrigToOrdered[i]);
205  }
206  }
207 
208 
209  // Walk sub-ordering (=column index) out.
210  labelList patchFaces(nFaces);
211  List<topoDistanceData> patchData(nFaces);
212  nFaces = 0;
213  forAllConstIter(labelHashSet, patchIDs, iter)
214  {
215  const polyPatch& pp = pbm[iter.key()];
216  const labelUList& fc = pp.faceCells();
217  forAll(fc, i)
218  {
219  patchFaces[nFaces] = pp.start()+i;
220  patchData[nFaces] = topoDistanceData
221  (
222  orderedToOld[fc[i]],// passive data: global column
223  0 // distance: layer
224  );
225  nFaces++;
226  }
227  }
228 
229  // Field on cells and faces.
230  List<topoDistanceData> cellData(mesh.nCells());
231  List<topoDistanceData> faceData(mesh.nFaces());
232 
233  // Propagate information inwards
235  (
236  mesh,
237  patchFaces,
238  patchData,
239  faceData,
240  cellData,
241  0
242  );
243 
244  deltaCalc.iterate(nLayers_);
245 
246  Info<< type() << " : did not visit "
247  << deltaCalc.getUnsetCells()
248  << " cells out of " << nTotalCells
249  << "; using " << method_().type() << " renumbering for these" << endl;
250 
251  // Get cell order using the method(). These values will get overwritten
252  // by any visited cell so are used only if the number of nLayers is limited.
253  labelList oldToOrdered
254  (
255  invert
256  (
257  mesh.nCells(),
258  method_().renumber(mesh, points)
259  )
260  );
261 
262  // Use specialised sorting to sorted either layers or columns first
263  // Done so that at no point we need to combine both into a single
264  // index and we might run out of label size.
266  (
267  cellData,
268  orderedToOld,
269  layerLess(depthFirst_, oldToOrdered, cellData)
270  );
271 
272  // Return furthest away cell first
273  if (reverse_)
274  {
275  reverse(orderedToOld);
276  }
277 
278  return orderedToOld;
279 }
280 
281 
282 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Macros for easy insertion into run-time selection tables.
label getUnsetCells() const
Get number of unvisited cells, i.e. cells that were not (yet)
Definition: FaceCellWave.C:778
virtual label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
Definition: FaceCellWave.C:951
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:109
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Version of FaceCellWave that walks through prismatic cells only.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Post-processing mesh subset tool. Given the original mesh and the list of selected cells,...
Definition: fvMeshSubset.H:74
const labelList & cellMap() const
Return cell map.
void setLargeCellSubset(const labelList &region, const label currentRegion, const label patchID=-1, const bool syncCouples=true)
Set the subset from all cells with region == currentRegion.
Definition: fvMeshSubset.C:893
const fvMesh & subMesh() const
Return reference to subset mesh.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:64
label toGlobal(const label i) const
From local to global.
Definition: globalIndexI.H:82
label nTotalCells() const
Return total number of cells in decomposed mesh.
Foam::polyBoundaryMesh.
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch indices corresponding to the given names.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1515
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:404
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:280
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:313
label nCells() const
label nFaces() const
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.....
structuredRenumber(const dictionary &renumberDict)
Construct given the renumber dictionary.
virtual labelList renumber(const pointField &) const
Return the order in which cells need to be visited, i.e.
For use with FaceCellWave. Determines topological distance to starting faces.
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
const pointField & points
volScalarField & b
Definition: createFields.H:25
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
labelList invert(const label len, const labelUList &)
Invert one-to-one map. Unmapped elements will be -1.
Definition: ListOps.C:37
ListType renumber(const labelUList &oldToNew, const ListType &)
Renumber the values (not the indices) of a list.
messageStream Info
void reverse(UList< T > &, const label n)
Definition: UListI.H:334
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
defineTypeNameAndDebug(combustionModel, 0)
error FatalError
static const label labelMax
Definition: label.H:62
void sortedOrder(const UList< T > &, labelList &order)
Generate the (stable) sort order for the list.
static const char nl
Definition: Ostream.H:266
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488