37Foam::activeBaffleVelocityFvPatchVectorField::
38activeBaffleVelocityFvPatchVectorField
44 fixedValueFvPatchVectorField(
p, iF),
47 cyclicPatchLabel_(-1),
54 maxOpenFractionDelta_(0),
59Foam::activeBaffleVelocityFvPatchVectorField::
60activeBaffleVelocityFvPatchVectorField
68 fixedValueFvPatchVectorField(ptf,
p, iF, mapper),
70 cyclicPatchName_(ptf.cyclicPatchName_),
71 cyclicPatchLabel_(ptf.cyclicPatchLabel_),
72 orientation_(ptf.orientation_),
73 initWallSf_(ptf.initWallSf_),
74 initCyclicSf_(ptf.initCyclicSf_),
75 nbrCyclicSf_(ptf.nbrCyclicSf_),
76 openFraction_(ptf.openFraction_),
77 openingTime_(ptf.openingTime_),
78 maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
83Foam::activeBaffleVelocityFvPatchVectorField::
84activeBaffleVelocityFvPatchVectorField
91 fixedValueFvPatchVectorField(
p, iF,
dict, false),
92 pName_(
dict.getOrDefault<
word>(
"p",
"p")),
94 cyclicPatchLabel_(
p.patch().
boundaryMesh().findPatchID(cyclicPatchName_)),
95 orientation_(
dict.get<label>(
"orientation")),
104 ).neighbFvPatch().Sf()
106 openFraction_(
dict.get<scalar>(
"openFraction")),
107 openingTime_(
dict.get<scalar>(
"openingTime")),
108 maxOpenFractionDelta_(
dict.get<scalar>(
"maxOpenFractionDelta")),
115Foam::activeBaffleVelocityFvPatchVectorField::
116activeBaffleVelocityFvPatchVectorField
121 fixedValueFvPatchVectorField(ptf),
123 cyclicPatchName_(ptf.cyclicPatchName_),
124 cyclicPatchLabel_(ptf.cyclicPatchLabel_),
125 orientation_(ptf.orientation_),
126 initWallSf_(ptf.initWallSf_),
127 initCyclicSf_(ptf.initCyclicSf_),
128 nbrCyclicSf_(ptf.nbrCyclicSf_),
129 openFraction_(ptf.openFraction_),
130 openingTime_(ptf.openingTime_),
131 maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
136Foam::activeBaffleVelocityFvPatchVectorField::
137activeBaffleVelocityFvPatchVectorField
143 fixedValueFvPatchVectorField(ptf, iF),
145 cyclicPatchName_(ptf.cyclicPatchName_),
146 cyclicPatchLabel_(ptf.cyclicPatchLabel_),
147 orientation_(ptf.orientation_),
148 initWallSf_(ptf.initWallSf_),
149 initCyclicSf_(ptf.initCyclicSf_),
150 nbrCyclicSf_(ptf.nbrCyclicSf_),
151 openFraction_(ptf.openFraction_),
152 openingTime_(ptf.openingTime_),
153 maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
165 fixedValueFvPatchVectorField::autoMap(m);
173 const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
174 initWallSf_ = patch().patchSlice(areas);
175 initCyclicSf_ = patch().boundaryMesh()
179 nbrCyclicSf_ = refCast<const cyclicFvPatch>
185 ).neighbFvPatch().patch().patchSlice(areas);
195 fixedValueFvPatchVectorField::rmap(ptf, addr);
198 const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
199 initWallSf_ = patch().patchSlice(areas);
200 initCyclicSf_ = patch().boundaryMesh()
204 nbrCyclicSf_ = refCast<const cyclicFvPatch>
210 ).neighbFvPatch().patch().patchSlice(areas);
222 if (curTimeIndex_ != this->db().time().
timeIndex())
229 const fvPatch& cyclicPatch = patch().boundaryMesh()[cyclicPatchLabel_];
230 const labelList& cyclicFaceCells = cyclicPatch.patch().faceCells();
231 const fvPatch& nbrPatch = refCast<const cyclicFvPatch>
237 scalar forceDiff = 0;
240 forAll(cyclicFaceCells, facei)
242 forceDiff +=
p[cyclicFaceCells[facei]]*
mag(initCyclicSf_[facei]);
246 forAll(nbrFaceCells, facei)
248 forceDiff -=
p[nbrFaceCells[facei]]*
mag(nbrCyclicSf_[facei]);
259 this->db().time().deltaTValue()/openingTime_,
260 maxOpenFractionDelta_
262 *(orientation_*
sign(forceDiff)),
268 Info<<
"openFraction = " << openFraction_ <<
endl;
271 const vectorField newSfw((1 - openFraction_)*initWallSf_);
274 Sfw[facei] = newSfw[facei];
280 openFraction_*initCyclicSf_;
282 mag(cyclicPatch.Sf());
285 openFraction_*nbrCyclicSf_;
289 curTimeIndex_ = this->db().time().timeIndex();
292 fixedValueFvPatchVectorField::updateCoeffs();
303 os.
writeEntry(
"maxOpenFractionDelta", maxOpenFractionDelta_);
305 writeEntry(
"value",
os);
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
SubField is a Field obtained as a section of another Field, without its own allocation....
This velocity boundary condition simulates the opening of a baffle due to local flow conditions,...
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void rmap(const fvPatchVectorField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
virtual void operator=(const UList< vector > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
const scalarField & magSf() const
Return face area magnitudes.
const polyPatch & patch() const
Return the polyPatch.
const vectorField & Sf() const
Return face area vectors.
const labelUList & faceCells() const
Return face-cell addressing.
Lookup type of boundary radiation properties.
A class for handling words, derived from Foam::string.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
To & refCast(From &r)
Reference type cast template function.
dimensionedScalar sign(const dimensionedScalar &ds)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
static constexpr const zero Zero
Global zero (0)
#define forAll(list, i)
Loop across all elements in list.