dynamicRefineFvMesh.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2017 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::dynamicRefineFvMesh
26 
27 Description
28  A fvMesh with built-in refinement.
29 
30  Determines which cells to refine/unrefine and does all in update().
31 
32 
33  // How often to refine
34  refineInterval 1;
35  // Field to be refinement on
36  field alpha.water;
37  // Refine field inbetween lower..upper
38  lowerRefineLevel 0.001;
39  upperRefineLevel 0.999;
40  // If value < unrefineLevel (default=GREAT) unrefine
41  //unrefineLevel 10;
42  // Have slower than 2:1 refinement
43  nBufferLayers 1;
44  // Refine cells only up to maxRefinement levels
45  maxRefinement 2;
46  // Stop refinement if maxCells reached
47  maxCells 200000;
48  // Flux field and corresponding velocity field. Fluxes on changed
49  // faces get recalculated by interpolating the velocity. Use 'none'
50  // on surfaceScalarFields that do not need to be reinterpolated, use
51  // NaN to detect use of mapped variable
52  correctFluxes
53  (
54  (phi none) //NaN) //none)
55  (nHatf none) //none)
56  (rho*phi none) //none)
57  (ghf none) //NaN) //none)
58  );
59  // Write the refinement level as a volScalarField
60  dumpLevel true;
61 
62 
63 SourceFiles
64  dynamicRefineFvMesh.C
65 
66 \*---------------------------------------------------------------------------*/
67 
68 #ifndef dynamicRefineFvMesh_H
69 #define dynamicRefineFvMesh_H
70 
71 #include "dynamicFvMesh.H"
72 #include "hexRef8.H"
73 #include "PackedBoolList.H"
74 #include "Switch.H"
75 
76 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
77 
78 namespace Foam
79 {
80 
81 /*---------------------------------------------------------------------------*\
82  Class dynamicRefineFvMesh Declaration
83 \*---------------------------------------------------------------------------*/
84 
86 :
87  public dynamicFvMesh
88 {
89 protected:
90 
91  //- Mesh cutting engine
93 
94  //- Dump cellLevel for postprocessing
96 
97  //- Fluxes to map
99 
100  //- Number of refinement/unrefinement steps done so far.
102 
103  //- Protected cells (usually since not hexes)
105 
106 
107  // Protected Member Functions
108 
109  //- Count set/unset elements in packedlist.
110  static label count(const PackedBoolList&, const unsigned int);
111 
112  //- Calculate cells that cannot be refined since would trigger
113  // refinement of protectedCell_ (since 2:1 refinement cascade)
114  void calculateProtectedCells(PackedBoolList& unrefineableCell) const;
115 
116  //- Read the projection parameters from dictionary
117  void readDict();
118 
119 
120  //- Refine cells. Update mesh and fields.
122 
123  //- Unrefine cells. Gets passed in centre points of cells to combine.
125 
126 
127  // Selection of cells to un/refine
128 
129  //- Calculates approximate value for refinement level so
130  // we don't go above maxCell
131  scalar getRefineLevel
132  (
133  const label maxCells,
134  const label maxRefinement,
135  const scalar refineLevel,
136  const scalarField&
137  ) const;
138 
139  //- Get per cell max of connected point
140  scalarField maxPointField(const scalarField&) const;
141 
142  //- Get point max of connected cell
144 
145  scalarField cellToPoint(const scalarField& vFld) const;
146 
148  (
149  const scalarField& fld,
150  const scalar minLevel,
151  const scalar maxLevel
152  ) const;
153 
154  //- Select candidate cells for refinement
155  virtual void selectRefineCandidates
156  (
157  const scalar lowerRefineLevel,
158  const scalar upperRefineLevel,
159  const scalarField& vFld,
160  PackedBoolList& candidateCell
161  ) const;
162 
163  //- Subset candidate cells for refinement
165  (
166  const label maxCells,
167  const label maxRefinement,
168  const PackedBoolList& candidateCell
169  ) const;
170 
171  //- Select points that can be unrefined.
173  (
174  const scalar unrefineLevel,
175  const PackedBoolList& markedCell,
176  const scalarField& pFld
177  ) const;
178 
179  //- Extend markedCell with cell-face-cell.
180  void extendMarkedCells(PackedBoolList& markedCell) const;
181 
182  //- Check all cells have 8 anchor points
184  (
186  label& nProtected
187  ) const;
188 
189 private:
190 
191  //- Disallow default bitwise copy construct
193 
194  //- Disallow default bitwise assignment
195  void operator=(const dynamicRefineFvMesh&);
196 
197 public:
198 
199  //- Runtime type information
200  TypeName("dynamicRefineFvMesh");
201 
202 
203  // Constructors
204 
205  //- Construct from IOobject
206  explicit dynamicRefineFvMesh(const IOobject& io);
207 
208 
209  //- Destructor
210  virtual ~dynamicRefineFvMesh();
211 
212 
213  // Member Functions
214 
215  //- Direct access to the refinement engine
216  const hexRef8& meshCutter() const
217  {
218  return meshCutter_;
219  }
220 
221  //- Cells which should not be refined/unrefined
222  const PackedBoolList& protectedCell() const
223  {
224  return protectedCell_;
225  }
226 
227  //- Cells which should not be refined/unrefined
229  {
230  return protectedCell_;
231  }
232 
233  //- Update the mesh for both mesh motion and topology change
234  virtual bool update();
235 
236 
237  // Writing
238 
239  //- Write using given format, version and compression
240  virtual bool writeObject
241  (
245  const bool valid
246  ) const;
247 
248 };
249 
250 
251 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
252 
253 } // End namespace Foam
254 
255 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
256 
257 #endif
258 
259 // ************************************************************************* //
virtual labelList selectRefineCells(const label maxCells, const label maxRefinement, const PackedBoolList &candidateCell) const
Subset candidate cells for refinement.
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
scalarField maxPointField(const scalarField &) const
Get per cell max of connected point.
static label count(const PackedBoolList &, const unsigned int)
Count set/unset elements in packedlist.
const hexRef8 & meshCutter() const
Direct access to the refinement engine.
PackedBoolList protectedCell_
Protected cells (usually since not hexes)
scalar getRefineLevel(const label maxCells, const label maxRefinement, const scalar refineLevel, const scalarField &) const
Calculates approximate value for refinement level so.
scalarField cellToPoint(const scalarField &vFld) const
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool valid) const
Write using given format, version and compression.
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
void checkEightAnchorPoints(PackedBoolList &protectedCell, label &nProtected) const
Check all cells have 8 anchor points.
HashTable< word > correctFluxes_
Fluxes to map.
TypeName("dynamicRefineFvMesh")
Runtime type information.
void calculateProtectedCells(PackedBoolList &unrefineableCell) const
Calculate cells that cannot be refined since would trigger.
A fvMesh with built-in refinement.
virtual bool update()
Update the mesh for both mesh motion and topology change.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
autoPtr< mapPolyMesh > unrefine(const labelList &)
Unrefine cells. Gets passed in centre points of cells to combine.
Refinement of (split) hexes using polyTopoChange.
Definition: hexRef8.H:64
hexRef8 meshCutter_
Mesh cutting engine.
void extendMarkedCells(PackedBoolList &markedCell) const
Extend markedCell with cell-face-cell.
const PackedBoolList & protectedCell() const
Cells which should not be refined/unrefined.
streamFormat
Enumeration for the format of data in the stream.
Definition: IOstream.H:86
autoPtr< mapPolyMesh > refine(const labelList &)
Refine cells. Update mesh and fields.
An STL-conforming hash table.
Definition: HashTable.H:62
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:193
label nRefinementIterations_
Number of refinement/unrefinement steps done so far.
virtual void selectRefineCandidates(const scalar lowerRefineLevel, const scalar upperRefineLevel, const scalarField &vFld, PackedBoolList &candidateCell) const
Select candidate cells for refinement.
void readDict()
Read the projection parameters from dictionary.
A bit-packed bool list.
virtual labelList selectUnrefinePoints(const scalar unrefineLevel, const PackedBoolList &markedCell, const scalarField &pFld) const
Select points that can be unrefined.
Abstract base class for geometry and/or topology changing fvMesh.
Definition: dynamicFvMesh.H:51
Version number type.
Definition: IOstream.H:96
virtual ~dynamicRefineFvMesh()
Destructor.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
scalarField error(const scalarField &fld, const scalar minLevel, const scalar maxLevel) const
Switch dumpLevel_
Dump cellLevel for postprocessing.
Namespace for OpenFOAM.
scalarField maxCellField(const volScalarField &) const
Get point max of connected cell.