38 const DimensionedField<scalar, volMesh>& iF
41 mixedFvPatchScalarField(p, iF),
44 psiName_(
"thermo:psi"),
46 accommodationCoeff_(1.0),
47 Twall_(p.size(), 0.0),
52 valueFraction() = 0.0;
59 const DimensionedField<scalar, volMesh>& iF,
60 const dictionary& dict
63 mixedFvPatchScalarField(p, iF),
64 UName_(dict.lookupOrDefault<word>(
"U",
"U")),
65 rhoName_(dict.lookupOrDefault<word>(
"rho",
"rho")),
66 psiName_(dict.lookupOrDefault<word>(
"psi",
"thermo:psi")),
67 muName_(dict.lookupOrDefault<word>(
"mu",
"thermo:mu")),
68 accommodationCoeff_(dict.
lookup<scalar>(
"accommodationCoeff")),
69 Twall_(
"Twall", dict, p.size()),
70 gamma_(dict.lookupOrDefault<scalar>(
"gamma", 1.4))
74 mag(accommodationCoeff_) < small
75 ||
mag(accommodationCoeff_) > 2.0
81 ) <<
"unphysical accommodationCoeff specified" 82 <<
"(0 < accommodationCoeff <= 1)" <<
endl 86 if (dict.found(
"value"))
100 valueFraction() = 0.0;
106 const smoluchowskiJumpTFvPatchScalarField& ptf,
108 const DimensionedField<scalar, volMesh>& iF,
109 const fvPatchFieldMapper& mapper
112 mixedFvPatchScalarField(ptf, p, iF, mapper),
114 rhoName_(ptf.rhoName_),
115 psiName_(ptf.psiName_),
116 muName_(ptf.muName_),
117 accommodationCoeff_(ptf.accommodationCoeff_),
118 Twall_(mapper(ptf.Twall_)),
125 const smoluchowskiJumpTFvPatchScalarField& ptpsf,
126 const DimensionedField<scalar, volMesh>& iF
129 mixedFvPatchScalarField(ptpsf, iF),
130 accommodationCoeff_(ptpsf.accommodationCoeff_),
131 Twall_(ptpsf.Twall_),
140 const fvPatchFieldMapper& m
143 mixedFvPatchScalarField::autoMap(m);
150 const fvPatchField<scalar>& ptf,
156 const smoluchowskiJumpTFvPatchScalarField& ptpsf =
157 refCast<const smoluchowskiJumpTFvPatchScalarField>(ptf);
159 Twall_.rmap(ptpsf.Twall_, addr);
174 const fvPatchField<scalar>& ppsi =
180 const dictionary& thermophysicalProperties =
187 thermophysicalProperties.subDict(
"mixture").subDict(
"transport")
195 *2.0*gamma_/Pr.
value()/(gamma_ + 1.0)
196 *(2.0 - accommodationCoeff_)/accommodationCoeff_
199 Field<scalar> aCoeff(prho.snGrad() - prho/C2);
200 Field<scalar> KEbyRho(0.5*
magSqr(pU));
202 valueFraction() = (1.0/(1.0 + patch().deltaCoeffs()*C2));
206 mixedFvPatchScalarField::updateCoeffs();
214 writeEntryIfDifferent<word>(os,
"U",
"U", UName_);
215 writeEntryIfDifferent<word>(os,
"rho",
"rho", rhoName_);
216 writeEntryIfDifferent<word>(os,
"psi",
"thermo:psi", psiName_);
217 writeEntryIfDifferent<word>(os,
"mu",
"thermo:mu", muName_);
219 writeEntry(os,
"accommodationCoeff", accommodationCoeff_);
233 smoluchowskiJumpTFvPatchScalarField
virtual void write(Ostream &) const
Write.
fvPatchField< vector > fvPatchVectorField
errorManipArg< error, int > exit(error &err, const int errNo=1)
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Vector< scalar > vector
A scalar version of the templated Vector.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void write(Ostream &) const
Write.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Macros for easy insertion into run-time selection tables.
const word dictName() const
Return the local dictionary name (final part of scoped name)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
stressControl lookup("compactNormalStress") >> compactNormalStress
fvPatchField< scalar > fvPatchScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const Type & value() const
Return const reference to value.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
List< label > labelList
A List of labels.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
smoluchowskiJumpTFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void operator=(const UList< Type > &)
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const scalar piByTwo(0.5 *pi)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)