37Foam::activePressureForceBaffleVelocityFvPatchVectorField::
38activePressureForceBaffleVelocityFvPatchVectorField
44 fixedValueFvPatchVectorField(
p, iF),
47 cyclicPatchLabel_(-1),
53 maxOpenFractionDelta_(0),
55 minThresholdValue_(0),
62Foam::activePressureForceBaffleVelocityFvPatchVectorField::
63activePressureForceBaffleVelocityFvPatchVectorField
71 fixedValueFvPatchVectorField(ptf,
p, iF, mapper),
73 cyclicPatchName_(ptf.cyclicPatchName_),
74 cyclicPatchLabel_(ptf.cyclicPatchLabel_),
75 initWallSf_(ptf.initWallSf_),
76 initCyclicSf_(ptf.initCyclicSf_),
77 nbrCyclicSf_(ptf.nbrCyclicSf_),
78 openFraction_(ptf.openFraction_),
79 openingTime_(ptf.openingTime_),
80 maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
82 minThresholdValue_(ptf.minThresholdValue_),
84 baffleActivated_(ptf.baffleActivated_),
85 opening_(ptf.opening_)
89Foam::activePressureForceBaffleVelocityFvPatchVectorField::
90activePressureForceBaffleVelocityFvPatchVectorField
97 fixedValueFvPatchVectorField(
p, iF,
dict, false),
98 pName_(
dict.getOrDefault<
word>(
"p",
"p")),
100 cyclicPatchLabel_(
p.patch().
boundaryMesh().findPatchID(cyclicPatchName_)),
104 openFraction_(
dict.get<scalar>(
"openFraction")),
105 openingTime_(
dict.get<scalar>(
"openingTime")),
106 maxOpenFractionDelta_(
dict.get<scalar>(
"maxOpenFractionDelta")),
108 minThresholdValue_(
dict.get<scalar>(
"minThresholdValue")),
109 fBased_(
dict.get<
bool>(
"forceBased")),
117 initWallSf_ =
p.Sf();
118 initCyclicSf_ =
p.boundaryMesh()[cyclicPatchLabel_].Sf();
119 nbrCyclicSf_ = refCast<const cyclicFvPatch>
121 p.boundaryMesh()[cyclicPatchLabel_],
123 ).neighbFvPatch().Sf();
130Foam::activePressureForceBaffleVelocityFvPatchVectorField::
131activePressureForceBaffleVelocityFvPatchVectorField
136 fixedValueFvPatchVectorField(ptf),
138 cyclicPatchName_(ptf.cyclicPatchName_),
139 cyclicPatchLabel_(ptf.cyclicPatchLabel_),
140 initWallSf_(ptf.initWallSf_),
141 initCyclicSf_(ptf.initCyclicSf_),
142 nbrCyclicSf_(ptf.nbrCyclicSf_),
143 openFraction_(ptf.openFraction_),
144 openingTime_(ptf.openingTime_),
145 maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
147 minThresholdValue_(ptf.minThresholdValue_),
148 fBased_(ptf.fBased_),
149 baffleActivated_(ptf.baffleActivated_),
150 opening_(ptf.opening_)
154Foam::activePressureForceBaffleVelocityFvPatchVectorField::
155activePressureForceBaffleVelocityFvPatchVectorField
161 fixedValueFvPatchVectorField(ptf, iF),
163 cyclicPatchName_(ptf.cyclicPatchName_),
164 cyclicPatchLabel_(ptf.cyclicPatchLabel_),
165 initWallSf_(ptf.initWallSf_),
166 initCyclicSf_(ptf.initCyclicSf_),
167 nbrCyclicSf_(ptf.nbrCyclicSf_),
168 openFraction_(ptf.openFraction_),
169 openingTime_(ptf.openingTime_),
170 maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
172 minThresholdValue_(ptf.minThresholdValue_),
173 fBased_(ptf.fBased_),
174 baffleActivated_(ptf.baffleActivated_),
175 opening_(ptf.opening_)
186 fixedValueFvPatchVectorField::autoMap(m);
198 Info <<
"faceArea[active] "<< i <<
endl;
202 if (patch().size() > 0)
204 const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
205 initWallSf_ = patch().patchSlice(areas);
207 patch().boundaryMesh()[cyclicPatchLabel_].patchSlice(areas);
208 nbrCyclicSf_ = refCast<const cyclicFvPatch>
214 ).neighbFvPatch().patch().patchSlice(areas);
225 fixedValueFvPatchVectorField::rmap(ptf, addr);
228 const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
229 initWallSf_ = patch().patchSlice(areas);
231 patch().boundaryMesh()[cyclicPatchLabel_].patchSlice(areas);
232 nbrCyclicSf_ = refCast<const cyclicFvPatch>
238 ).neighbFvPatch().patch().patchSlice(areas);
250 if (curTimeIndex_ != this->db().time().
timeIndex())
255 const fvPatch& cyclicPatch = patch().boundaryMesh()[cyclicPatchLabel_];
256 const labelList& cyclicFaceCells = cyclicPatch.patch().faceCells();
258 refCast<const cyclicFvPatch>(cyclicPatch).neighbFvPatch();
262 scalar valueDiff = 0;
266 forAll(cyclicFaceCells, facei)
268 valueDiff +=
p[cyclicFaceCells[facei]]*
mag(initCyclicSf_[facei]);
269 area +=
mag(initCyclicSf_[facei]);
273 forAll(nbrFaceCells, facei)
275 valueDiff -=
p[nbrFaceCells[facei]]*
mag(initCyclicSf_[facei]);
280 valueDiff = valueDiff/(area + VSMALL);
289 Info<<
"Force difference (threshold) = " << valueDiff
290 <<
"(" << minThresholdValue_ <<
")" <<
endl;
294 Info<<
"Area-averaged pressure difference (threshold) = "
295 << valueDiff <<
"(" << minThresholdValue_ <<
")" <<
endl;
299 if (
mag(valueDiff) >
mag(minThresholdValue_) || baffleActivated_)
309 this->db().time().deltaT().value()/openingTime_,
310 maxOpenFractionDelta_
317 baffleActivated_ =
true;
321 openFraction_ =
max(
min(1 - 1
e-6, openFraction_), 1
e-6);
326 Info<<
"Open fraction = " << openFraction_ <<
endl;
329 scalar areaFraction = 0.0;
333 areaFraction = openFraction_;
337 areaFraction = 1 - openFraction_;
346 initWallSf_ = patch().Sf();
347 initCyclicSf_ = patch().boundaryMesh()[cyclicPatchLabel_].Sf();
348 nbrCyclicSf_ = refCast<const cyclicFvPatch>
351 ).neighbFvPatch().Sf();
356 vectorField newSfw((1 - areaFraction)*initWallSf_);
359 Sfw[facei] = newSfw[facei];
364 const_cast<vectorField&
>(cyclicPatch.Sf()) = areaFraction*initCyclicSf_;
366 const_cast<scalarField&
>(cyclicPatch.magSf()) =
mag(cyclicPatch.Sf());
369 const_cast<vectorField&
>(nbrPatch.
Sf()) = areaFraction*nbrCyclicSf_;
373 curTimeIndex_ = this->db().time().timeIndex();
376 fixedValueFvPatchVectorField::updateCoeffs();
387 os.
writeEntry(
"maxOpenFractionDelta", maxOpenFractionDelta_);
392 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....
void size(const label n)
Older name for setAddressableSize.
This boundary condition is applied to the flow velocity, to simulate the opening or closure of a baff...
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,...
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
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.
splitCell * master() const
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.
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)
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
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.