45 void Foam::mixerFvMesh::addZonesAndModifiers()
54 || topoChanger_.size()
58 <<
"Zones and modifiers already present. Skipping." 64 Info<<
"Time = " << time().timeName() << endl
65 <<
"Adding zones and modifiers to the mesh" <<
endl;
68 List<pointZone*> pz(1);
83 List<faceZone*> fz(3);
86 const word innerSliderName(motionDict_.subDict(
"slider").lookup(
"inside"));
87 const polyPatch& innerSlider = boundaryMesh()[innerSliderName];
93 isf[i] = innerSlider.start() + i;
100 boolList(innerSlider.size(),
false),
106 const word outerSliderName(motionDict_.subDict(
"slider").lookup(
"outside"));
107 const polyPatch& outerSlider = boundaryMesh()[outerSliderName];
113 osf[i] = outerSlider.start() + i;
120 boolList(outerSlider.size(),
false),
135 List<cellZone*> cz(1);
138 regionSplit rs(*
this);
141 label originRegion = rs[findNearestCell(cs().origin())];
144 label nMovingCells = 0;
148 if (rs[celli] == originRegion)
150 movingCells[nMovingCells] = celli;
155 movingCells.setSize(nMovingCells);
156 Info<<
"Number of cells in the moving region: " << nMovingCells <<
endl;
166 Info<<
"Adding point, face and cell zones" <<
endl;
167 addZones(pz, fz, cz);
170 Info<<
"Adding topology modifiers" <<
endl;
171 topoChanger_.setSize(1);
180 outerSliderName +
"Zone",
181 innerSliderName +
"Zone",
195 void Foam::mixerFvMesh::calcMovingMasks()
const 200 <<
"Calculating point and cell masks" 204 if (movingPointsMaskPtr_)
207 <<
"point mask already calculated" 213 scalarField& movingPointsMask = *movingPointsMaskPtr_;
218 const labelList& cellAddr = cellZones()[
"movingCells"];
222 const cell& curCell = c[cellAddr[celli]];
227 const face& curFace = f[curCell[facei]];
231 movingPointsMask[curFace[pointi]] = 1;
236 const word innerSliderZoneName
238 word(motionDict_.subDict(
"slider").lookup(
"inside"))
242 const labelList& innerSliderAddr = faceZones()[innerSliderZoneName];
244 forAll(innerSliderAddr, facei)
246 const face& curFace = f[innerSliderAddr[facei]];
250 movingPointsMask[curFace[pointi]] = 1;
254 const word outerSliderZoneName
256 word(motionDict_.subDict(
"slider").lookup(
"outside"))
260 const labelList& outerSliderAddr = faceZones()[outerSliderZoneName];
262 forAll(outerSliderAddr, facei)
264 const face& curFace = f[outerSliderAddr[facei]];
268 movingPointsMask[curFace[pointi]] = 0;
277 Foam::mixerFvMesh::mixerFvMesh
296 ).optionalSubDict(typeName +
"Coeffs")
303 motionDict_.subDict(
"coordinateSystem")
307 movingPointsMaskPtr_(
nullptr)
309 addZonesAndModifiers();
311 Info<<
"Mixer mesh:" <<
nl 312 <<
" origin: " << cs().origin() <<
nl 313 <<
" axis: " << cs().axis() <<
nl 314 <<
" rpm: " << rpm_ <<
endl;
330 if (!movingPointsMaskPtr_)
335 return *movingPointsMaskPtr_;
344 csPtr_->globalPosition
346 csPtr_->localPosition(
points())
347 +
vector(0, rpm_*360.0*time().deltaTValue()/60.0, 0)
355 if (topoChangeMap.
valid())
367 csPtr_->globalPosition
369 csPtr_->localPosition(oldPoints())
370 +
vector(0, rpm_*360.0*time().deltaTValue()/60.0, 0)
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Vector< scalar > vector
A scalar version of the templated Vector.
Macros for easy insertion into run-time selection tables.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
List< bool > boolList
Bool container classes.
virtual ~mixerFvMesh()
Destructor.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Abstract base class for a topology changing fvMesh.
List< label > labelList
A List of labels.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
bool valid() const
Return true if the autoPtr valid (ie, the pointer is set)
errorManip< error > abort(error &err)
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
static autoPtr< coordinateSystem > New(const objectRegistry &obr, const dictionary &dict)
Select constructed from dictionary and objectRegistry.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
virtual bool update()
Update the mesh for both mesh motion and topology change.
void deleteDemandDrivenData(DataPtr &dataPtr)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
List< cell > cellList
list of cells
#define InfoInFunction
Report an information message using Foam::Info.