Go to the documentation of this file.
53 const Foam::scalar Foam::layerAdditionRemoval::addDelta_ = 0.3;
54 const Foam::scalar Foam::layerAdditionRemoval::removeDelta_ = 0.1;
59 void Foam::layerAdditionRemoval::checkDefinition()
64 <<
"Master face zone named " << faceZoneID_.
name()
65 <<
" cannot be found."
71 minLayerThickness_ < VSMALL
72 || maxLayerThickness_ < minLayerThickness_
76 <<
"Incorrect layer thickness definition."
82 reduce(nFaces, sumOp<label>());
87 <<
"Face extrusion zone contains no faces. "
88 <<
"Please check your mesh definition."
94 Pout<<
"Cell layer addition/removal object " <<
name() <<
" :" <<
nl
95 <<
" faceZoneID: " << faceZoneID_ <<
endl;
100 Foam::scalar Foam::layerAdditionRemoval::readOldThickness
102 const dictionary&
dict
109 void Foam::layerAdditionRemoval::clearAddressing()
const
111 if (pointsPairingPtr_)
115 Pout<<
"layerAdditionRemoval::clearAddressing()" <<
nl
116 <<
" clearing pointsPairingPtr_" <<
endl;
122 if (facesPairingPtr_)
126 Pout<<
"layerAdditionRemoval::clearAddressing()" <<
nl
127 <<
" clearing facesPairingPtr_" <<
endl;
137 Foam::layerAdditionRemoval::layerAdditionRemoval
142 const word& zoneName,
143 const scalar minThickness,
144 const scalar maxThickness,
145 const bool thicknessFromVolume
150 minLayerThickness_(minThickness),
151 maxLayerThickness_(maxThickness),
152 thicknessFromVolume_(thicknessFromVolume),
153 oldLayerThickness_(-1.0),
154 pointsPairingPtr_(
nullptr),
155 facesPairingPtr_(
nullptr),
163 Foam::layerAdditionRemoval::layerAdditionRemoval
173 minLayerThickness_(
dict.
get<scalar>(
"minLayerThickness")),
174 maxLayerThickness_(
dict.
get<scalar>(
"maxLayerThickness")),
176 oldLayerThickness_(readOldThickness(
dict)),
177 pointsPairingPtr_(
nullptr),
178 facesPairingPtr_(
nullptr),
199 if (triggerRemoval_ > -1 || triggerAddition_ > -1)
220 if (
min(V) < -VSMALL)
223 <<
"negative cell volume. Error in mesh motion before "
224 <<
"topological change.\n V: " << V
229 scalar minDelta = GREAT;
233 if (thicknessFromVolume_)
238 scalar curDelta = V[mc[facei]]/
mag(S[fz[facei]]);
239 avgDelta += curDelta;
240 minDelta =
min(minDelta, curDelta);
241 maxDelta =
max(maxDelta, curDelta);
249 const Map<label>& zoneMeshPointMap = fz().meshPointMap();
259 const edge&
e = cellEdges[i];
261 if (zoneMeshPointMap.found(
e[0]))
263 if (!zoneMeshPointMap.found(
e[1]))
266 avgDelta += curDelta;
268 minDelta =
min(minDelta, curDelta);
269 maxDelta =
max(maxDelta, curDelta);
274 if (zoneMeshPointMap.found(
e[1]))
277 avgDelta += curDelta;
279 minDelta =
min(minDelta, curDelta);
280 maxDelta =
max(maxDelta, curDelta);
296 Pout<<
"bool layerAdditionRemoval::changeTopology() const "
297 <<
" for object " <<
name() <<
" : " <<
nl
298 <<
"Layer thickness: min: " << minDelta
299 <<
" max: " << maxDelta <<
" avg: " << avgDelta
300 <<
" old thickness: " << oldLayerThickness_ <<
nl
301 <<
"Removal threshold: " << minLayerThickness_
302 <<
" addition threshold: " << maxLayerThickness_ <<
endl;
305 bool topologicalChange =
false;
309 if (oldLayerThickness_ < 0)
313 Pout<<
"First step. No addition/removal" <<
endl;
317 oldLayerThickness_ = avgDelta;
319 topologicalChange =
false;
321 else if (avgDelta < oldLayerThickness_)
324 if (minDelta < minLayerThickness_)
327 if (setLayerPairing())
339 Pout<<
"bool layerAdditionRemoval::changeTopology() "
340 <<
" const for object " <<
name() <<
" : "
341 <<
"Triggering layer removal" <<
endl;
348 oldLayerThickness_ = GREAT;
350 topologicalChange =
true;
361 oldLayerThickness_ = avgDelta;
367 if (maxDelta > maxLayerThickness_)
371 Pout<<
"bool layerAdditionRemoval::changeTopology() const "
372 <<
" for object " <<
name() <<
" : "
373 <<
"Triggering layer addition" <<
endl;
380 oldLayerThickness_ = 0;
382 topologicalChange =
true;
386 oldLayerThickness_ = avgDelta;
390 return topologicalChange;
399 if (triggerRemoval_ == topoChanger().
mesh().time().
timeIndex())
401 removeCellLayer(
ref);
406 Pout<<
"layerAdditionRemoval::setRefinement(polyTopoChange&) "
407 <<
"for object " <<
name() <<
" : "
408 <<
"Clearing addressing after layer removal" <<
endl;
411 triggerRemoval_ = -1;
415 if (triggerAddition_ == topoChanger().
mesh().time().
timeIndex())
422 Pout<<
"layerAdditionRemoval::setRefinement(polyTopoChange&) "
423 <<
"for object " <<
name() <<
" : "
424 <<
"Clearing addressing after layer addition" <<
endl;
427 triggerAddition_ = -1;
437 Pout<<
"layerAdditionRemoval::updateMesh(const mapPolyMesh&) "
438 <<
"for object " <<
name() <<
" : "
439 <<
"Clearing addressing on external request";
441 if (pointsPairingPtr_ || facesPairingPtr_)
452 faceZoneID_.update(topoChanger().
mesh().faceZones());
460 if (t < VSMALL || maxLayerThickness_ < t)
463 <<
"Incorrect layer thickness definition."
467 minLayerThickness_ = t;
473 if (t < minLayerThickness_)
476 <<
"Incorrect layer thickness definition."
480 maxLayerThickness_ = t;
489 << minLayerThickness_ <<
nl
490 << oldLayerThickness_ <<
nl
491 << maxLayerThickness_ <<
nl
492 << thicknessFromVolume_ <<
endl;
499 <<
" type " <<
type()
501 <<
" faceZoneName " << faceZoneID_.name()
503 <<
" minLayerThickness " << minLayerThickness_
505 <<
" maxLayerThickness " << maxLayerThickness_
507 <<
" thicknessFromVolume " << thicknessFromVolume_
509 <<
" oldLayerThickness " << oldLayerThickness_
511 <<
" active " << active()
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
virtual const pointField & points() const
Return raw points.
virtual bool changeTopology() const
Check for topology change.
A class for handling words, derived from Foam::string.
virtual void updateMesh(const mapPolyMesh &)
Force recalculation of locally stored data on topological change.
const polyTopoChanger & topoChanger() const
Return reference to morph engine.
List of mesh modifiers defining the mesh dynamics.
const cellList & cells() const
Direct mesh changes based on v1.3 polyTopoChange syntax.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
virtual void setRefinement(polyTopoChange &) const
Insert the layer addition/removal instructions.
Ostream & endl(Ostream &os)
Add newline and flush stream.
edgeList edges(const faceUList &f) const
Return cell edges.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
T lookupOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
const faceZoneMesh & faceZones() const
Return face zone mesh.
virtual ~layerAdditionRemoval()
Destructor.
void deleteDemandDrivenData(DataPtr &dataPtr)
void setMaxLayerThickness(const scalar t) const
Set max layer thickness which triggers removal.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
A subset of mesh faces organised as a primitive patch.
word name(const complex &c)
Return string representation of complex.
void setMinLayerThickness(const scalar t) const
Set min layer thickness which triggers removal.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
ITstream & lookup(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
errorManip< error > abort(error &err)
const scalarField & cellVolumes() const
Begin block [isseparator].
Virtual base class for mesh modifiers.
const word & name() const
Return name of this modifier.
const labelList & masterCells() const
virtual const faceList & faces() const
Return raw faces.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const polyMesh & mesh() const
Return the mesh reference.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const dimensionedScalar e
Elementary charge.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
virtual void write(Ostream &) const
Write.
const Time & time() const
Return the top-level database.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
label timeIndex() const
Return current time index.
const keyType & name() const
Return name.
virtual void writeDict(Ostream &) const
Write dictionary.
defineTypeNameAndDebug(combustionModel, 0)
A cell is defined as a list of faces with extra functionality.
bool active() const
Has the zone been found.
label index() const
Return index of first matching zone.
const vectorField & faceAreas() const