38Foam::outletMachNumberPressureFvPatchScalarField::
39outletMachNumberPressureFvPatchScalarField
45 fixedValueFvPatchScalarField(
p, iF),
58Foam::outletMachNumberPressureFvPatchScalarField::
59outletMachNumberPressureFvPatchScalarField
66 fixedValueFvPatchScalarField(
p, iF,
dict),
67 M_(
dict.getOrDefault<scalar>(
"M", 0)),
68 pBack_(
dict.get<scalar>(
"pBack")),
69 c1_(
dict.getOrDefault<scalar>(
"c1", 0)),
70 A1_(
dict.getOrDefault<scalar>(
"A1", 0)),
71 phiName_(
dict.getOrDefault<
word>(
"phi",
"phi")),
72 rhoName_(
dict.getOrDefault<
word>(
"rho",
"rho")),
73 UName_(
dict.getOrDefault<
word>(
"U",
"U")),
75 relax_(
dict.getOrDefault<scalar>(
"relax", 0))
79Foam::outletMachNumberPressureFvPatchScalarField::
80outletMachNumberPressureFvPatchScalarField
88 fixedValueFvPatchScalarField(ptf,
p, iF, mapper),
93 phiName_(ptf.phiName_),
94 rhoName_(ptf.rhoName_),
101Foam::outletMachNumberPressureFvPatchScalarField::
102outletMachNumberPressureFvPatchScalarField
107 fixedValueFvPatchScalarField(tppsf),
109 pBack_(tppsf.pBack_),
112 phiName_(tppsf.phiName_),
113 rhoName_(tppsf.rhoName_),
114 UName_(tppsf.UName_),
115 choked_(tppsf.choked_),
120Foam::outletMachNumberPressureFvPatchScalarField::
121outletMachNumberPressureFvPatchScalarField
127 fixedValueFvPatchScalarField(tppsf, iF),
129 pBack_(tppsf.pBack_),
132 phiName_(tppsf.phiName_),
133 rhoName_(tppsf.rhoName_),
134 UName_(tppsf.UName_),
135 choked_(tppsf.choked_),
152 this->internalField().name()
155 const label patchi = patch().index();
175 const vectorField UbOld(
U.oldTime().boundaryField()[patchi]);
178 Ub = relax_*UbOld + (1 - relax_)*Ub;
203 <<
"Pelase specify M in the dictionary"
210 if (A1_ == 0.0 || c1_ == 0.0)
217 const scalar area =
gSum(
mag(patch().Sf()));
224 if (
M[i] < 0 || r[i] >=1)
227 <<
"or pBack/ptot ratio is larger then one"
239 pbNew = relax_*pb + (1 -relax_)*pbNew;
243 fixedValueFvPatchScalarField::updateCoeffs();
261 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...
const dimensionSet & dimensions() const
Return dimensions.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const Internal & internalField() const
Return a const-reference to the dimensioned internal field.
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.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
static const word dictName
virtual tmp< volScalarField > gamma() const =0
Gamma = Cp/Cv [].
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Fundamental fluid thermodynamic properties.
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
This boundary condition maintains a certain subsonic Mach number at an outlet patch by dynamically ad...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A class for handling words, derived from Foam::string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
#define WarningInFunction
Report a warning using Foam::Warning.
Type gSum(const FieldField< Field, Type > &f)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const dimensionSet dimArea(sqr(dimLength))
const dimensionSet dimVelocity
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Type gAverage(const FieldField< Field, Type > &f)
const dimensionSet dimDensity
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define forAll(list, i)
Loop across all elements in list.