62bool limitRefinementLevel
71 label oldNCells = refCells.
size();
75 const labelList& cCells = cellCells[celli];
79 if (refLevel[cCells[i]] > (refLevel[celli]+1))
89 if (refCells.
size() > oldNCells)
91 Info<<
"Added an additional " << refCells.
size() - oldNCells
92 <<
" cells to satisfy 1:2 refinement level"
102int main(
int argc,
char *argv[])
106 "Attempt to determine refinement levels of a refined cartesian mesh.\n"
107 "Run BEFORE snapping!"
110 argList::addBoolOption
113 "Read level from refinementLevel file"
120 Info<<
"Dividing cells into bins depending on cell volume.\nThis will"
121 <<
" correspond to refinement levels for a mesh with only 2x2x2"
123 <<
"The upper range for every bin is always 1.1 times the lower range"
124 <<
" to allow for some truncation error."
127 const bool readLevel =
args.
found(
"readLevel");
142 lowerLimits.
append(sortedVols[0]);
147 if (sortedVols[i] > upperLimits.
last())
156 Info<<
"Collected " << bin.
size() <<
" elements in bin "
157 << lowerLimits.
last() <<
" .. "
162 lowerLimits.
append(sortedVols[i]);
165 Info<<
"Creating new bin " << lowerLimits.
last()
166 <<
" .. " << upperLimits.
last()
173 bin.
append(sortedVols.indices()[i]);
177 bins.
last().shrink();
187 Info<<
"Volume bins:" <<
nl;
195 Info<<
" " << lowerLimits[binI] <<
" .. " << upperLimits[binI]
196 <<
" : writing " << bin.
size() <<
" cells to cellSet "
215 polyMesh::defaultRegion,
233 p[patchi] =
patches[patchi].clone(fMesh.boundaryMesh()).ptr();
236 fMesh.addFvPatches(
p);
244 polyMesh::defaultRegion,
248 if (!readLevel && refHeader.typeHeaderOk<
labelIOList>(
true))
251 <<
"Detected " << refHeader.name() <<
" file in "
252 << polyMesh::defaultRegion <<
" directory. Please remove to"
253 <<
" recreate it or use the -readLevel option to use it"
299 refLevel[bin[i]] = bins.
size() - binI - 1;
300 postRefLevel[bin[i]] = refLevel[bin[i]];
305 postRefLevel.boundaryFieldRef();
310 forAll(postRefLevel.boundaryField(), patchi)
316 Info<<
"Setting field for patch "<<
endl;
320 label own =
mesh.faceOwner()[pp.
start() + facei];
322 bField[facei] = postRefLevel[own];
326 Info<<
"Determined current refinement level and writing to "
327 << postRefLevel.name() <<
" (as volScalarField; for post processing)"
329 << polyMesh::defaultRegion/refLevel.name()
330 <<
" (as labelIOList; for meshing)" <<
nl
334 postRefLevel.write();
355 Info<<
"Collected " << refCells.
size() <<
" cells that need to be"
356 <<
" refined to get closer to overall 2:1 refinement level limit"
358 <<
"Written cells to be refined to cellSet " << refCells.
name()
363 Info<<
"After refinement this tool can be run again to see if the 2:1"
364 <<
" limit is observed all over the mesh" <<
nl <<
endl;
368 Info<<
"All cells in the mesh observe the 2:1 refinement level limit"
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
DynamicList< T, SizeMin > & shrink()
Shrink the allocated space to the number of elements used.
void append(const T &val)
Copy append an element to the end of this list.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
label size() const noexcept
The number of elements in table.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const word & name() const noexcept
Return the object name.
A list that is sorted upon construction or when explicitly requested with the sort() method.
void size(const label n)
Older name for setAddressableSize.
T & last()
Return the last element of the list.
bool found(const word &optName) const
Return true if the named option is found.
A collection of cell labels.
Mesh data needed to do the Finite Volume discretisation.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
A patch is a list of labels that address the faces in the global face list.
label start() const
Return start label of this patch in the polyMesh face list.
Cell-face mesh analysis engine.
virtual bool write(const bool valid=true) const
Write using setting from DB.
const polyBoundaryMesh & patches
#define WarningInFunction
Report a warning using Foam::Warning.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
constexpr char nl
The newline '\n' character (0x0a)
Foam::argList args(argc, argv)
#define forAll(list, i)
Loop across all elements in list.