61 volumetric_(
dict.getOrDefault(
"volumetric", true))
69 rhoName_ =
dict.getOrDefault<word>(
"rho",
"rho");
73 if (!this->readValueEntry(
dict))
83 const matchedFlowRateOutletVelocityFvPatchVectorField& ptf,
85 const DimensionedField<vector, volMesh>& iF,
86 const fvPatchFieldMapper& mapper
89 fixedValueFvPatchField<
vector>(ptf,
p, iF, mapper),
90 inletPatchName_(ptf.inletPatchName_),
91 rhoName_(ptf.rhoName_),
92 volumetric_(ptf.volumetric_)
103 inletPatchName_(ptf.inletPatchName_),
104 rhoName_(ptf.rhoName_),
105 volumetric_(ptf.volumetric_)
117 inletPatchName_(ptf.inletPatchName_),
118 rhoName_(ptf.rhoName_),
119 volumetric_(ptf.volumetric_)
125 template<
class RhoType>
126 void Foam::matchedFlowRateOutletVelocityFvPatchVectorField::updateValues
128 const label inletPatchID,
129 const RhoType& rhoOutlet,
130 const RhoType& rhoInlet
133 const fvPatch&
p =
patch();
134 const fvPatch& inletPatch =
p.boundaryMesh()[inletPatchID];
148 nUp =
max(nUp, scalar(0));
153 const_cast<volVectorField&>
155 dynamic_cast<const volVectorField&>(internalField())
163 inletPatchU.updateCoeffs();
166 const scalar flowRate = -
gSum(rhoInlet*(inletPatch.Sf() & inletPatchU));
169 const scalar estimatedFlowRate =
gSum(rhoOutlet*(
patch().magSf()*nUp));
171 if (estimatedFlowRate > 0.5*flowRate)
173 nUp *= (
mag(flowRate)/
mag(estimatedFlowRate));
177 nUp += ((flowRate - estimatedFlowRate)/
gSum(rhoOutlet*
patch().magSf()));
196 const label inletPatchID =
197 patch().patch().boundaryMesh().findPatchID(inletPatchName_);
199 if (inletPatchID < 0)
202 <<
"Unable to find inlet patch " << inletPatchName_
208 updateValues(inletPatchID, one{}, one{});
213 if (db().foundObject<volScalarField>(rhoName_))
223 rho.boundaryField()[
patch().index()],
224 rho.boundaryField()[inletPatchID]
234 fixedValueFvPatchVectorField::updateCoeffs();
261 matchedFlowRateOutletVelocityFvPatchVectorField
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
fvPatchField< vector > fvPatchVectorField
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
matchedFlowRateOutletVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
Velocity outlet boundary condition which corrects the extrapolated velocity to match the flow rate of...
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
void writeValueEntry(Ostream &os) const
Write *this field as a "value" entry.
virtual void write(Ostream &) const
Write.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Type gSum(const FieldField< Field, Type > &f)
virtual void write(Ostream &) const
Write.
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
string evaluate(label fieldWidth, const std::string &s, size_t pos=0, size_t len=std::string::npos)
String evaluation with specified (positive, non-zero) field width.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
OBJstream os(runTime.globalPath()/outputName)
A simple container of IOobject preferences. Can also be used for general handling of read/no-read/rea...
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const std::string patch
OpenFOAM patch number as a std::string.
Field< vector > vectorField
Specialisation of Field<T> for vector.
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
"buffered" : (MPI_Bsend, MPI_Recv)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)