48 void Foam::linearValveFvMesh::addZonesAndModifiers()
61 <<
"Zones and modifiers already present. Skipping." 68 <<
"Adding zones and modifiers to the mesh" <<
endl;
71 List<pointZone*> pz(1);
74 pz[0] =
new pointZone(
"cutPointZone", 0,
pointZones());
79 List<faceZone*> fz(3);
82 const word innerSliderName
84 motionDict_.
subDict(
"slider").
get<word>(
"inside")
86 const polyPatch& innerSlider =
boundaryMesh()[innerSliderName];
98 const word outerSliderName
100 motionDict_.
subDict(
"slider").
get<word>(
"outside")
102 const polyPatch& outerSlider =
boundaryMesh()[outerSliderName];
114 fz[2] =
new faceZone(
"cutFaceZone", 2,
faceZones());
116 List<cellZone*> cz(0);
118 Info<<
"Adding point, face and cell zones" <<
endl;
122 Info<<
"Adding topology modifiers" <<
endl;
132 outerSliderName +
"Zone",
133 innerSliderName +
"Zone",
149 void Foam::linearValveFvMesh::makeSlidersDead()
151 const polyTopoChanger& topoChanges = topoChanger_;
156 if (isA<slidingInterface>(topoChanges[modI]))
158 topoChanges[modI].disable();
163 <<
"Don't know what to do with mesh modifier " 164 << modI <<
" of type " << topoChanges[modI].type()
171 void Foam::linearValveFvMesh::makeSlidersLive()
173 const polyTopoChanger& topoChanges = topoChanger_;
178 if (isA<slidingInterface>(topoChanges[modI]))
180 topoChanges[modI].enable();
185 <<
"Don't know what to do with mesh modifier " 186 << modI <<
" of type " << topoChanges[modI].type()
193 bool Foam::linearValveFvMesh::attached()
const 195 const polyTopoChanger& topoChanges = topoChanger_;
201 if (isA<slidingInterface>(topoChanges[modI]))
205 || refCast<const slidingInterface>(topoChanges[modI]).attached();
212 if (isA<slidingInterface>(topoChanges[modI]))
217 != refCast<const slidingInterface>(topoChanges[modI]).attached()
222 <<
" named " << topoChanges[modI].name()
223 <<
" out of sync: Should be" << result
231 Info<<
"linearValveFvMesh: attached!" <<
endl;
235 Info<<
"linearValveFvMesh: detached!" <<
endl;
245 Foam::linearValveFvMesh::linearValveFvMesh(
const IOobject&
io)
259 ).optionalSubDict(typeName +
"Coeffs")
263 addZonesAndModifiers();
280 Info<<
"Decoupling sliding interfaces" <<
endl;
288 msPtr_->updateMesh();
292 Info<<
"Sliding interfaces decoupled" <<
endl;
300 setMorphTimeIndex(3*time().
timeIndex() + 1);
303 msPtr_->updateMesh();
307 if (topoChangeMap().hasMotionPoints())
309 Info<<
"Topology change; executing pre-motion" <<
endl;
310 movePoints(topoChangeMap().preMotionPoints());
317 movePoints(msPtr_->curPoints());
320 Info<<
"Coupling sliding interfaces" <<
endl;
323 setMorphTimeIndex(3*time().
timeIndex() + 2);
326 Info<<
"Moving points post slider attach" <<
endl;
328 msPtr_->updateMesh();
330 Info<<
"Sliding interfaces coupled: " << attached() <<
endl;
virtual bool update()
Update the mesh for both mesh motion and topology change.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
Virtual base class for mesh motion solver.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T. FatalIOError if not found, or if the number of tokens is incorrect.
const Time & time() const
Return the top-level database.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Macros for easy insertion into run-time selection tables.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
#define forAll(list, i)
Loop across all elements in list.
label size() const noexcept
The number of elements in table.
writeOption writeOpt() const noexcept
Get the write option.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
label size() const noexcept
The number of entries in the list.
void addZones(PtrList< pointZone > &&pz, PtrList< faceZone > &&fz, PtrList< cellZone > &&cz)
Add mesh zones.
Abstract base class for a topology changing fvMesh.
errorManip< error > abort(error &err)
const T * set(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
virtual bool write(const bool writeOnProc=true) const
Write mesh using IO settings from time.
static word timeName(const scalar t, const int precision=precision_)
Return time name of given scalar time formatted with the given precision.
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
defineTypeNameAndDebug(combustionModel, 0)
void setSize(const label newLen)
Same as resize()
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
const pointZoneMesh & pointZones() const noexcept
Return point zone mesh.
virtual ~linearValveFvMesh()
Destructor.
Automatically write from objectRegistry::writeObject()
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
messageStream Info
Information stream (stdout output on master, null elsewhere)
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, IOobject::NO_REGISTER)
polyTopoChanger topoChanger_
Defines the attributes of an object for which implicit objectRegistry management is supported...
#define InfoInFunction
Report an information message using Foam::Info.