fvMeshTopoChangersRefiner.H
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) 2011-2023 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 Class
25  Foam::fvMeshTopoChangers::refiner
26 
27 Description
28  Dynamic mesh refinement/unrefinement based on volScalarField values.
29 
30  Refinement can optionally be specified in a cellZone or in multiple
31  regions, each controlled by a different volScalarField.
32 
33 Usage
34  Example of single field based refinement in all cells:
35  \verbatim
36  topoChanger
37  {
38  type refiner;
39 
40  libs ("libfvMeshTopoChangers.so");
41 
42  // How often to refine
43  refineInterval 1;
44 
45  // Field to be refinement on
46  field alpha.water;
47 
48  // Refine field in between lower..upper
49  lowerRefineLevel 0.001;
50  upperRefineLevel 0.999;
51 
52  // Have slower than 2:1 refinement
53  nBufferLayers 1;
54 
55  // Refine cells only up to maxRefinement levels
56  maxRefinement 2;
57 
58  // Stop refinement if maxCells reached
59  maxCells 200000;
60 
61  // Flux field and corresponding velocity field. Fluxes on changed
62  // faces get recalculated by interpolating the velocity. Use 'none'
63  // on surfaceScalarFields that do not need to be reinterpolated.
64  correctFluxes
65  (
66  (phi none)
67  (nHatf none)
68  (rhoPhi none)
69  (alphaPhi.water none)
70  (ghf none)
71  );
72 
73  // Write the refinement level as a volScalarField
74  dumpLevel true;
75  }
76  \endverbatim
77 
78  Example of single field based refinement in two regions:
79  \verbatim
80  topoChanger
81  {
82  type refiner;
83 
84  // How often to refine
85  refineInterval 1;
86 
87  refinementRegions
88  {
89  region1
90  {
91  cellZone refinementRegion1;
92 
93  // Field to be refinement on
94  field alpha.water;
95 
96  // Refine field in between lower..upper
97  lowerRefineLevel 0.001;
98  upperRefineLevel 0.999;
99 
100  // Refine cells only up to maxRefinement levels
101  maxRefinement 1;
102 
103  // If value < unrefineLevel unrefine
104  unrefineLevel 10;
105  }
106 
107  region2
108  {
109  cellZone refinementRegion2;
110 
111  // Field to be refinement on
112  field alpha.water;
113 
114  // Refine field in between lower..upper
115  lowerRefineLevel 0.001;
116  upperRefineLevel 0.999;
117 
118  // Refine cells only up to maxRefinement levels
119  maxRefinement 2;
120 
121  // If value < unrefineLevel unrefine
122  unrefineLevel 10;
123  }
124  }
125 
126  // If value < unrefineLevel (default=great) unrefine
127  // unrefineLevel 10;
128  // Have slower than 2:1 refinement
129  nBufferLayers 1;
130 
131  // Stop refinement if maxCells reached
132  maxCells 200000;
133 
134  // Flux field and corresponding velocity field. Fluxes on changed
135  // faces get recalculated by interpolating the velocity. Use 'none'
136  // on surfaceScalarFields that do not need to be reinterpolated.
137  correctFluxes
138  (
139  (phi none)
140  (nHatf none)
141  (rhoPhi none)
142  (alphaPhi.water none)
143  (ghf none)
144  );
145 
146  // Write the refinement level as a volScalarField
147  dumpLevel true;
148  }
149  \endverbatim
150 
151 
152 SourceFiles
153  fvMeshTopoChangersRefiner.C
154 
155 \*---------------------------------------------------------------------------*/
156 
157 #ifndef fvMeshTopoChangersRefiner_H
158 #define fvMeshTopoChangersRefiner_H
159 
160 #include "fvMeshTopoChanger.H"
161 #include "hexRef8.H"
162 #include "PackedBoolList.H"
163 #include "Switch.H"
164 
165 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
166 
167 namespace Foam
168 {
169 namespace fvMeshTopoChangers
170 {
171 
172 /*---------------------------------------------------------------------------*\
173  Class fvMeshTopoChangers::refiner Declaration
174 \*---------------------------------------------------------------------------*/
175 
176 class refiner
177 :
178  public fvMeshTopoChanger
179 {
180  // Private Member Data
181 
182  //- Refinement control dictionary
183  dictionary dict_;
184 
185  label refineInterval_;
186 
187  label maxCells_;
188 
189  label nBufferLayers_;
190 
191  //- Mesh cutting engine
192  hexRef8 meshCutter_;
193 
194  //- Dump cellLevel for postprocessing
195  Switch dumpLevel_;
196 
197  //- Fluxes to map
198  HashTable<word> correctFluxes_;
199 
200  //- Number of refinement/unrefinement steps done so far.
201  label nRefinementIterations_;
202 
203  //- Protected cells (usually since not hexes)
204  PackedBoolList protectedCells_;
205 
206  //- Has the mesh refinement changed since the last time
207  // the refinement history was written
208  mutable bool changedSinceWrite_;
209 
210  //- The time index used for updating
211  label timeIndex_;
212 
213 
214  // Private Member Functions
215 
216  //- Count set/unset elements in packedlist.
217  static label count(const PackedBoolList&, const unsigned int);
218 
219  //- Calculate cells that cannot be refined since would trigger
220  // refinement of protectedCells_ (since 2:1 refinement cascade)
221  void calculateProtectedCells(PackedBoolList& unrefineableCells) const;
222 
223  //- Read the parameters from dictionary
224  void readDict();
225 
226  //- Refine cells. Update mesh and fields.
227  autoPtr<polyTopoChangeMap> refine(const labelList&);
228 
229  //- Unrefine cells. Gets passed in centre points of cells to combine.
230  autoPtr<polyTopoChangeMap> unrefine(const labelList&);
231 
232  //- Find the U field name corresponding to Uf
233  word Uname(const surfaceVectorField& Uf) const;
234 
235  void refineFluxes
236  (
237  const labelHashSet& masterFaces,
238  const polyTopoChangeMap& map
239  );
240 
241  void unrefineFluxes
242  (
243  const Map<label>& faceToSplitPoint,
244  const polyTopoChangeMap& map
245  );
246 
247  void refineUfs
248  (
249  const labelHashSet& masterFaces,
250  const polyTopoChangeMap& map
251  );
252 
253  void unrefineUfs
254  (
255  const Map<label>& faceToSplitPoint,
256  const polyTopoChangeMap& map
257  );
258 
259 
260  // Selection of cells to un/refine
261 
262  const cellZone& findCellZone
263  (
264  const word& cellZoneName
265  ) const;
266 
267  scalarField cellToPoint(const scalarField& vFld) const;
268 
270  (
271  const scalarField& fld,
272  const scalar minLevel,
273  const scalar maxLevel
274  ) const;
275 
277  (
278  const scalarField& fld,
279  const labelList& cells,
280  const scalar minLevel,
281  const scalar maxLevel
282  ) const;
283 
284  //- Select candidate cells for refinement
285  virtual void selectRefineCandidates
286  (
287  PackedBoolList& candidateCell,
288  const scalar lowerRefineLevel,
289  const scalar upperRefineLevel,
290  const scalar maxRefinement,
291  const scalarField& vFld
292  ) const;
293 
294  //- Select candidate cells for refinement
295  virtual void selectRefineCandidates
296  (
297  PackedBoolList& candidateCell,
298  const scalar lowerRefineLevel,
299  const scalar upperRefineLevel,
300  const scalar maxRefinement,
301  const scalarField& vFld,
302  const labelList& cells
303  ) const;
304 
305  //- Select candidate cells for refinement
306  virtual scalar selectRefineCandidates
307  (
308  PackedBoolList& candidateCell,
309  const dictionary& refineDict
310  ) const;
311 
312  //- Subset candidate cells for refinement
313  virtual labelList selectRefineCells
314  (
315  const label maxCells,
316  const label maxRefinement,
317  const PackedBoolList& candidateCell
318  ) const;
319 
320  //- Select points that can be unrefined.
321  virtual labelList selectUnrefinePoints
322  (
323  const PackedBoolList& markedCell
324  ) const;
325 
326  //- Extend markedCell with cell-face-cell.
327  void extendMarkedCells(PackedBoolList& markedCell) const;
328 
329  //- Check all cells have 8 anchor points
330  void checkEightAnchorPoints
331  (
333  label& nProtected
334  ) const;
335 
336 
337 public:
338 
339  //- Runtime type information
340  TypeName("refiner");
341 
342 
343  // Constructors
344 
345  //- Construct from fvMesh and dictionary
346  refiner(fvMesh& mesh, const dictionary& dict);
347 
348  //- Disallow default bitwise copy construction
349  refiner(const refiner&) = delete;
350 
351 
352  //- Destructor
353  virtual ~refiner();
354 
355 
356  // Member Functions
357 
358  //- Direct access to the refinement engine
359  const hexRef8& meshCutter() const
360  {
361  return meshCutter_;
362  }
363 
364  //- Cells which should not be refined/unrefined
365  const PackedBoolList& protectedCell() const
366  {
367  return protectedCells_;
368  }
369 
370  //- Cells which should not be refined/unrefined
372  {
373  return protectedCells_;
374  }
375 
376  //- Update the mesh for both mesh motion and topology change
377  virtual bool update();
378 
379  //- Update corresponding to the given map
380  virtual void topoChange(const polyTopoChangeMap&);
381 
382  //- Update from another mesh using the given map
383  virtual void mapMesh(const polyMeshMap&);
384 
385  //- Update corresponding to the given distribution map
386  virtual void distribute(const polyDistributionMap&);
387 
388 
389  // Writing
390 
391  //- Write using given format, version and compression
392  virtual bool write(const bool write = true) const;
393 
394 
395  // Member Operators
396 
397  //- Disallow default bitwise assignment
398  void operator=(const refiner&) = delete;
399 };
400 
401 
402 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
403 
404 } // End namespace fvMeshTopoChangers
405 } // End namespace Foam
406 
407 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
408 
409 #endif
410 
411 // ************************************************************************* //
Generic GeometricField class.
An STL-conforming hash table.
Definition: HashTable.H:127
A bit-packed bool list.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A topoSetSource to select points based on usage in cells.
Definition: cellToPoint.H:52
A subset of mesh cells.
Definition: cellZone.H:58
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Class to handle errors and exceptions in a simple, consistent stream-based manner.
Definition: error.H:71
Abstract base class for fvMesh topology changers.
fvMesh & mesh()
Return the fvMesh.
Dynamic mesh refinement/unrefinement based on volScalarField values.
void operator=(const refiner &)=delete
Disallow default bitwise assignment.
virtual void topoChange(const polyTopoChangeMap &)
Update corresponding to the given map.
TypeName("refiner")
Runtime type information.
virtual void distribute(const polyDistributionMap &)
Update corresponding to the given distribution map.
const PackedBoolList & protectedCell() const
Cells which should not be refined/unrefined.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
virtual bool write(const bool write=true) const
Write using given format, version and compression.
virtual bool update()
Update the mesh for both mesh motion and topology change.
const hexRef8 & meshCutter() const
Direct access to the refinement engine.
refiner(fvMesh &mesh, const dictionary &dict)
Construct from fvMesh and dictionary.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Refinement of (split) hexes using polyTopoChange.
Definition: hexRef8.H:65
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
A class for handling words, derived from string.
Definition: word.H:62
const cellShapeList & cells
gmvFile<< "tracers "<< particles.size()<< nl;{ pointField positions(particles.size());label particlei=0;forAllConstIter(Cloud< passiveParticle >, particles, iter) { positions[particlei++]=iter().position(mesh);} for(i=0;i< pTraits< point >::nComponents;i++) { forAll(positions, particlei) { gmvFile<< component(positions[particlei], i)<< ' ';} gmvFile<< nl;}}forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Namespace for OpenFOAM.
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
dictionary dict