Go to the documentation of this file.
49 void Foam::mixerFvMesh::addZonesAndModifiers()
62 <<
"Zones and modifiers already present. Skipping."
69 <<
"Adding zones and modifiers to the mesh" <<
endl;
72 List<pointZone*> pz(1);
75 pz[0] =
new pointZone(
"cutPointZone", 0,
pointZones());
80 List<faceZone*> fz(3);
83 const word innerSliderName
85 motionDict_.
subDict(
"slider").
get<word>(
"inside")
87 const polyPatch& innerSlider =
boundaryMesh()[innerSliderName];
92 identity(innerSlider.size(), innerSlider.start()),
99 const word outerSliderName
101 motionDict_.
subDict(
"slider").
get<word>(
"outside")
103 const polyPatch& outerSlider =
boundaryMesh()[outerSliderName];
108 identity(outerSlider.size(), outerSlider.start()),
115 fz[2] =
new faceZone(
"cutFaceZone", 2,
faceZones());
117 List<cellZone*> cz(1);
120 regionSplit rs(*
this);
126 label nMovingCells = 0;
130 if (rs[celli] == originRegion)
132 movingCells[nMovingCells] = celli;
137 movingCells.resize(nMovingCells);
138 Info<<
"Number of cells in the moving region: " << nMovingCells <<
endl;
143 std::move(movingCells),
148 Info<<
"Adding point, face and cell zones" <<
endl;
152 Info<<
"Adding topology modifiers" <<
endl;
162 outerSliderName +
"Zone",
163 innerSliderName +
"Zone",
177 void Foam::mixerFvMesh::calcMovingMasks()
const
181 if (movingPointsMaskPtr_)
184 <<
"point mask already calculated"
190 scalarField& movingPointsMask = *movingPointsMaskPtr_;
195 const labelList& cellAddr = cellZones()[
"movingCells"];
197 for (
const label celli : cellAddr)
199 const cell& curCell =
c[celli];
201 for (
const label facei : curCell)
204 const face& curFace =
f[facei];
208 movingPointsMask[curFace[pointi]] = 1;
213 const word innerSliderZoneName
215 motionDict_.subDict(
"slider").get<word>(
"inside") +
"Zone"
218 const labelList& innerSliderAddr = faceZones()[innerSliderZoneName];
220 for (
const label facei : innerSliderAddr)
222 const face& curFace =
f[facei];
226 movingPointsMask[curFace[pointi]] = 1;
230 const word outerSliderZoneName
232 motionDict_.subDict(
"slider").get<word>(
"outside") +
"Zone"
235 const labelList& outerSliderAddr = faceZones()[outerSliderZoneName];
237 for (
const label facei : outerSliderAddr)
239 const face& curFace =
f[facei];
243 movingPointsMask[curFace[pointi]] = 0;
251 Foam::mixerFvMesh::mixerFvMesh
270 ).optionalSubDict(typeName +
"Coeffs")
273 rpm_(motionDict_.get<scalar>(
"rpm")),
274 movingPointsMaskPtr_(
nullptr)
276 if (motionDict_.found(coordinateSystem::typeName_()))
287 addZonesAndModifiers();
289 Info<<
"Mixer mesh:" <<
nl
290 <<
" origin: " << csys_.origin() <<
nl
291 <<
" axis: " << csys_.e3() <<
nl
292 <<
" rpm: " << rpm_ <<
endl;
309 if (!movingPointsMaskPtr_)
314 return *movingPointsMaskPtr_;
327 csys_.localPosition(
points())
336 if (topoChangeMap.
valid())
347 csys_.localPosition(oldPoints())
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
List< label > labelList
A List of labels.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
label size() const noexcept
The number of elements in table.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
#define InfoInFunction
Report an information message using Foam::Info.
virtual bool write(const bool valid=true) const
Write mesh using IO settings from time.
static constexpr const zero Zero
Global zero (0)
Template functions to aid in the implementation of demand driven data.
static word timeName(const scalar t, const int precision=precision_)
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Unit conversion functions.
bool valid() const noexcept
True if the managed pointer is non-null.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
#define forAll(list, i)
Loop across all elements in list.
const faceZoneMesh & faceZones() const
Return face zone mesh.
void deleteDemandDrivenData(DataPtr &dataPtr)
virtual const point & origin() const
Return origin.
label findNearestCell(const point &location) const
Find the cell with the nearest cell centre to location.
label nCells() const
Number of mesh cells.
writeOption writeOpt() const
The write option.
messageStream Info
Information stream (uses stdout - output is on the master only)
#define DebugInFunction
Report an information message using Foam::Info.
const pointZoneMesh & pointZones() const
Return point zone mesh.
static autoPtr< coordinateSystem > New(word modelType, const objectRegistry &obr, const dictionary &dict)
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
List< cell > cellList
A List of cells.
void setSize(const label newLen)
Same as resize()
A cylindrical coordinate system (r-theta-z). The coordinate system angle theta is always in radians.
Macros for easy insertion into run-time selection tables.
errorManip< error > abort(error &err)
const T * set(const label i) const
Return const pointer to element (if set) or nullptr.
virtual bool update()
Update the mesh for both mesh motion and topology change.
constexpr scalar rpmToRads(const scalar rpm) noexcept
Conversion from revolutions/minute to radians/sec.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
List< face > faceList
A List of faces.
Abstract base class for a topology changing fvMesh.
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
const dimensionedScalar c
Speed of light in a vacuum.
const Time & time() const
Return the top-level database.
virtual ~mixerFvMesh()
Destructor.
constant condensation/saturation model.
defineTypeNameAndDebug(combustionModel, 0)
void addZones(const List< pointZone * > &pz, const List< faceZone * > &fz, const List< cellZone * > &cz)
Add mesh zones.
Base class for coordinate system specification, the default coordinate system type is cartesian .
polyTopoChanger topoChanger_