41 const DimensionedField<scalar, volMesh>& iF
44 mixedFvPatchScalarField(
p, iF),
47 psiName_(
"thermo:psi"),
49 accommodationCoeff_(1.0),
50 Twall_(
p.size(),
Zero),
55 valueFraction() = 0.0;
61 const smoluchowskiJumpTFvPatchScalarField& ptf,
63 const DimensionedField<scalar, volMesh>& iF,
64 const fvPatchFieldMapper& mapper
67 mixedFvPatchScalarField(ptf,
p, iF, mapper),
69 rhoName_(ptf.rhoName_),
70 psiName_(ptf.psiName_),
72 accommodationCoeff_(ptf.accommodationCoeff_),
81 const DimensionedField<scalar, volMesh>& iF,
82 const dictionary&
dict 85 mixedFvPatchScalarField(
p, iF),
86 UName_(
dict.getOrDefault<word>(
"U",
"U")),
87 rhoName_(
dict.getOrDefault<word>(
"rho",
"rho")),
88 psiName_(
dict.getOrDefault<word>(
"psi",
"thermo:psi")),
89 muName_(
dict.getOrDefault<word>(
"mu",
"thermo:mu")),
90 accommodationCoeff_(
dict.
get<scalar>(
"accommodationCoeff")),
91 Twall_(
"Twall",
dict,
p.size()),
92 gamma_(
dict.getOrDefault<scalar>(
"gamma", 1.4))
96 mag(accommodationCoeff_) < SMALL
97 ||
mag(accommodationCoeff_) > 2.0
101 <<
"unphysical accommodationCoeff specified" 102 <<
"(0 < accommodationCoeff <= 1)" <<
endl 106 if (!this->readValueEntry(
dict))
114 valueFraction() = 0.0;
120 const smoluchowskiJumpTFvPatchScalarField& ptpsf,
121 const DimensionedField<scalar, volMesh>& iF
124 mixedFvPatchScalarField(ptpsf, iF),
125 accommodationCoeff_(ptpsf.accommodationCoeff_),
126 Twall_(ptpsf.Twall_),
136 const fvPatchFieldMapper& m
139 mixedFvPatchScalarField::autoMap(m);
146 const fvPatchField<scalar>& ptf,
175 thermophysicalProperties.subDict(
"mixture").subDict(
"transport")
182 *2.0*gamma_/
Pr.
value()/(gamma_ + 1.0)
183 *(2.0 - accommodationCoeff_)/accommodationCoeff_
186 Field<scalar> aCoeff(prho.snGrad() - prho/C2);
187 Field<scalar> KEbyRho(0.5*
magSqr(pU));
189 valueFraction() = (1.0/(1.0 +
patch().deltaCoeffs()*C2));
193 mixedFvPatchScalarField::updateCoeffs();
207 os.
writeEntry(
"accommodationCoeff", accommodationCoeff_);
221 smoluchowskiJumpTFvPatchScalarField
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
word dictName() const
The local dictionary name (final part of scoped name)
dimensionedScalar Pr("Pr", dimless, laminarTransport)
const Type & value() const noexcept
Return const reference to value.
virtual void write(Ostream &) const
Write.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
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.
const dimensionSet dimless
Dimensionless.
virtual void write(Ostream &) const
Write.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Macros for easy insertion into run-time selection tables.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
fvPatchField< scalar > fvPatchScalarField
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
constexpr scalar piByTwo(0.5 *M_PI)
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch.
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
OBJstream os(runTime.globalPath()/outputName)
smoluchowskiJumpTFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const std::string patch
OpenFOAM patch number as a std::string.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
List< label > labelList
A List of labels.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...
static constexpr const zero Zero
Global zero (0)