Go to the documentation of this file.
48 void Foam::MRFZone::setMRFFaces()
68 if (cellZoneID_ != -1)
73 zoneCell[cellLabels[i]] =
true;
82 if (zoneCell[own[facei]] || zoneCell[nei[facei]])
94 const polyPatch& pp =
patches[patchi];
96 if (pp.coupled() || excludedPatches.found(patchi))
100 label facei = pp.start()+i;
102 if (zoneCell[own[facei]])
109 else if (!isA<emptyPolyPatch>(pp))
113 label facei = pp.start()+i;
115 if (zoneCell[own[facei]])
139 if (faceType[facei] == 1)
141 internalFaces_[nInternal++] = facei;
144 internalFaces_.
setSize(nInternal);
151 const polyPatch& pp =
patches[patchi];
155 label facei = pp.start() + patchFacei;
157 if (faceType[facei] == 1)
159 nIncludedFaces[patchi]++;
161 else if (faceType[facei] == 2)
163 nExcludedFaces[patchi]++;
170 forAll(nIncludedFaces, patchi)
172 includedFaces_[patchi].
setSize(nIncludedFaces[patchi]);
173 excludedFaces_[patchi].
setSize(nExcludedFaces[patchi]);
180 const polyPatch& pp =
patches[patchi];
184 label facei = pp.start() + patchFacei;
186 if (faceType[facei] == 1)
188 includedFaces_[patchi][nIncludedFaces[patchi]++] = patchFacei;
190 else if (faceType[facei] == 2)
192 excludedFaces_[patchi][nExcludedFaces[patchi]++] = patchFacei;
200 faceSet internalFaces(mesh_,
"internalFaces", internalFaces_);
201 Pout<<
"Writing " << internalFaces.size()
202 <<
" internal faces in MRF zone to faceSet "
204 internalFaces.
write();
206 faceSet MRFFaces(mesh_,
"includedFaces", 100);
207 forAll(includedFaces_, patchi)
209 forAll(includedFaces_[patchi], i)
211 label patchFacei = includedFaces_[patchi][i];
215 Pout<<
"Writing " << MRFFaces.size()
216 <<
" patch faces in MRF zone to faceSet "
220 faceSet excludedFaces(mesh_,
"excludedFaces", 100);
221 forAll(excludedFaces_, patchi)
223 forAll(excludedFaces_[patchi], i)
225 label patchFacei = excludedFaces_[patchi][i];
226 excludedFaces.insert(
patches[patchi].
start()+patchFacei);
229 Pout<<
"Writing " << excludedFaces.size()
230 <<
" faces in MRF zone with special handling to faceSet "
232 excludedFaces.
write();
239 Foam::MRFZone::MRFZone
244 const word& cellZoneName
250 active_(coeffs_.lookupOrDefault(
"active",
true)),
251 cellZoneName_(cellZoneName),
257 origin_(coeffs_.lookup(
"origin")),
258 axis_(coeffs_.lookup(
"axis")),
263 coeffs_.readEntry(
"cellZone", cellZoneName_);
272 cellZoneID_ = mesh_.cellZones().findZoneID(cellZoneName_);
274 axis_ = axis_/
mag(axis_);
278 mesh_.boundaryMesh().patchSet(excludedPatchNames_)
281 excludedPatchLabels_.setSize(excludedPatchSet.size());
284 for (
const label patchi : excludedPatchSet)
286 excludedPatchLabels_[i++] = patchi;
289 bool cellZoneFound = (cellZoneID_ != -1);
296 <<
"cannot find MRF cellZone " << cellZoneName_
309 return omega_->value(mesh_.time().timeOutputValue())*axis_;
319 if (cellZoneID_ == -1)
328 const vector Omega = this->Omega();
333 ddtUc[celli] += (Omega ^ Uc[celli]);
340 if (cellZoneID_ == -1)
350 const vector Omega = this->Omega();
357 Usource[celli] += V[celli]*(Omega ^
U[celli]);
365 Usource[celli] -= V[celli]*(Omega ^
U[celli]);
378 if (cellZoneID_ == -1)
388 const vector Omega = this->Omega();
395 Usource[celli] += V[celli]*
rho[celli]*(Omega ^
U[celli]);
403 Usource[celli] -= V[celli]*
rho[celli]*(Omega ^
U[celli]);
411 if (cellZoneID_ == -1)
418 const vector Omega = this->Omega();
425 U[celli] -= (Omega ^ (
C[celli] - origin_));
430 volVectorField::Boundary& Ubf =
U.boundaryFieldRef();
432 forAll(includedFaces_, patchi)
434 forAll(includedFaces_[patchi], i)
436 label patchFacei = includedFaces_[patchi][i];
437 Ubf[patchi][patchFacei] =
Zero;
442 forAll(excludedFaces_, patchi)
444 forAll(excludedFaces_[patchi], i)
446 label patchFacei = excludedFaces_[patchi][i];
447 Ubf[patchi][patchFacei] -=
449 ^ (
C.boundaryField()[patchi][patchFacei] - origin_));
479 makeRelativeRhoFlux(
rho,
phi);
485 if (cellZoneID_ == -1)
492 const vector Omega = this->Omega();
499 U[celli] += (Omega ^ (
C[celli] - origin_));
503 volVectorField::Boundary& Ubf =
U.boundaryFieldRef();
505 forAll(includedFaces_, patchi)
507 forAll(includedFaces_[patchi], i)
509 label patchFacei = includedFaces_[patchi][i];
510 Ubf[patchi][patchFacei] =
511 (Omega ^ (
C.boundaryField()[patchi][patchFacei] - origin_));
516 forAll(excludedFaces_, patchi)
518 forAll(excludedFaces_[patchi], i)
520 label patchFacei = excludedFaces_[patchi][i];
521 Ubf[patchi][patchFacei] +=
522 (Omega ^ (
C.boundaryField()[patchi][patchFacei] - origin_));
540 makeAbsoluteRhoFlux(
rho,
phi);
551 const vector Omega = this->Omega();
554 volVectorField::Boundary& Ubf =
U.boundaryFieldRef();
556 forAll(includedFaces_, patchi)
558 const vectorField& patchC = mesh_.Cf().boundaryField()[patchi];
562 forAll(includedFaces_[patchi], i)
564 label patchFacei = includedFaces_[patchi][i];
566 pfld[patchFacei] = (Omega ^ (patchC[patchFacei] - origin_));
583 omega_->writeData(os);
585 if (excludedPatchNames_.size())
587 os.
writeEntry(
"nonRotatingPatches", excludedPatchNames_);
599 coeffs_.readEntry(
"cellZone", cellZoneName_);
600 cellZoneID_ = mesh_.cellZones().findZoneID(cellZoneName_);
608 if (mesh_.topoChanging())
int debug
Static debugging option.
List< label > labelList
A List of labels.
A class representing the concept of a field of oneFields used to avoid unnecessary manipulations for ...
A class for handling words, derived from Foam::string.
A field of fields is a PtrList of fields with reference counting.
static constexpr const zero Zero
Global zero.
label nInternalFaces() const
Number of internal faces.
label nFaces() const
Number of mesh faces.
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
void makeAbsolute(volVectorField &U) const
Make the given relative velocity absolute within the MRF region.
const cellZoneMesh & cellZones() const
Return cell zone mesh.
List< bool > boolList
A List of bools.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual Ostream & beginBlock(const keyType &kw)
Write begin block group with the given name.
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.
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
#define forAll(list, i)
Loop across all elements in list.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
bool read(const dictionary &dict)
Read MRF dictionary.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
label nCells() const
Number of mesh cells.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
word name(const complex &c)
Return string representation of complex.
virtual const fileName & name() const
Return the name of the stream.
virtual const labelList & faceOwner() const
Return face owner.
virtual Ostream & endBlock()
Write end block group.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
Mesh data needed to do the Finite Volume discretisation.
void writeData(Ostream &os) const
Write.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
vector Omega() const
Return the current Omega vector.
void correctBoundaryVelocity(volVectorField &U) const
Correct the boundary velocity for the rotation of the MRF region.
void makeRelative(volVectorField &U) const
Make the given absolute velocity relative within the MRF region.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label ListType::const_reference const label start
A class representing the concept of a field of 1 used to avoid unnecessary manipulations for objects ...
A List of wordRe with additional matching capabilities.
static const word null
An empty word.
const polyBoundaryMesh & patches
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Graphite solid properties.
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys and label hasher.
void setSize(const label newSize)
Alias for resize(const label)
defineTypeNameAndDebug(combustionModel, 0)
void update()
Update MRFZone faces if the mesh topology changes.
virtual const labelList & faceNeighbour() const
Return face neighbour.
void addCoriolis(const volVectorField &U, volVectorField &ddtU) const
Add the Coriolis force contribution to the acceleration field.