FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 SourceFiles
33  dynamicRefineFvMesh.C
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #ifndef dynamicRefineFvMesh_H
38 #define dynamicRefineFvMesh_H
39 
41 #include <dynamicMesh/hexRef8.H>
43 #include <OpenFOAM/Switch.H>
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 namespace Foam
48 {
49 
50 /*---------------------------------------------------------------------------*\
51  Class dynamicRefineFvMesh Declaration
52 \*---------------------------------------------------------------------------*/
53 
55 :
56  public dynamicFvMesh
57 {
58 protected:
59 
60  //- Mesh cutting engine
62 
63  //- Dump cellLevel for postprocessing
65 
66  //- Fluxes to map
68 
69  //- Number of refinement/unrefinement steps done so far.
71 
72  //- Protected cells (usually since not hexes)
74 
75 
76  // Private Member Functions
77 
78  //- Count set/unset elements in packedlist.
79  static label count(const PackedBoolList&, const unsigned int);
80 
81  //- Calculate cells that cannot be refined since would trigger
82  // refinement of protectedCell_ (since 2:1 refinement cascade)
83  void calculateProtectedCells(PackedBoolList& unrefineableCell) const;
84 
85  //- Read the projection parameters from dictionary
86  void readDict();
87 
88 
89  //- Refine cells. Update mesh and fields.
91 
92  //- Unrefine cells. Gets passed in centre points of cells to combine.
94 
95 
96  // Selection of cells to un/refine
97 
98  //- Calculates approximate value for refinement level so
99  // we don't go above maxCell
100  scalar getRefineLevel
101  (
102  const label maxCells,
103  const label maxRefinement,
104  const scalar refineLevel,
105  const scalarField&
106  ) const;
107 
108  //- Get per cell max of connected point
109  scalarField maxPointField(const scalarField&) const;
110 
111  //- Get point min of connected cell
113 
114  scalarField cellToPoint(const scalarField& vFld) const;
115 
117  (
118  const scalarField& fld,
119  const scalar minLevel,
120  const scalar maxLevel
121  ) const;
122 
123  //- Select candidate cells for refinement
124  virtual void selectRefineCandidates
125  (
126  const scalar lowerRefineLevel,
127  const scalar upperRefineLevel,
128  const scalarField& vFld,
129  PackedBoolList& candidateCell
130  ) const;
131 
132  //- Subset candidate cells for refinement
134  (
135  const label maxCells,
136  const label maxRefinement,
137  const PackedBoolList& candidateCell
138  ) const;
139 
140  //- Select points that can be unrefined.
142  (
143  const scalar unrefineLevel,
144  const PackedBoolList& markedCell,
145  const scalarField& pFld
146  ) const;
147 
148  //- Extend markedCell with cell-face-cell.
149  void extendMarkedCells(PackedBoolList& markedCell) const;
150 
151 
152 private:
153 
154  //- Disallow default bitwise copy construct
156 
157  //- Disallow default bitwise assignment
158  void operator=(const dynamicRefineFvMesh&);
159 
160 public:
161 
162  //- Runtime type information
163  TypeName("dynamicRefineFvMesh");
164 
165 
166  // Constructors
167 
168  //- Construct from IOobject
169  explicit dynamicRefineFvMesh(const IOobject& io);
170 
171 
172  // Destructor
173 
174  virtual ~dynamicRefineFvMesh();
175 
176 
177  // Member Functions
178 
179  //- Direct access to the refinement engine
180  const hexRef8& meshCutter() const
181  {
182  return meshCutter_;
183  }
184 
185  //- Cells which should not be refined/unrefined
187  {
188  return protectedCell_;
189  }
190 
191  //- Cells which should not be refined/unrefined
193  {
194  return protectedCell_;
195  }
196 
197  //- Update the mesh for both mesh motion and topology change
198  virtual bool update();
199 
200 
201  // Writing
202 
203  //- Write using given format, version and compression
204  virtual bool writeObject
205  (
209  ) const;
210 
211 };
212 
213 
214 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
215 
216 } // End namespace Foam
217 
218 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
219 
220 #endif
221 
222 // ************************ vim: set sw=4 sts=4 et: ************************ //