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
103 <<
exit(FatalIOError);
106 if (
dict.found(
"value"))
108 fvPatchField<scalar>::operator=
115 fvPatchField<scalar>::operator=(patchInternalField());
120 valueFraction() = 0.0;
126 const smoluchowskiJumpTFvPatchScalarField& ptpsf,
127 const DimensionedField<scalar, volMesh>& iF
130 mixedFvPatchScalarField(ptpsf, iF),
131 accommodationCoeff_(ptpsf.accommodationCoeff_),
132 Twall_(ptpsf.Twall_),
142 const fvPatchFieldMapper& m
145 mixedFvPatchScalarField::autoMap(m);
152 const fvPatchField<scalar>& ptf,
172 const fvPatchField<scalar>& ppsi =
178 const dictionary& thermophysicalProperties =
185 thermophysicalProperties.subDict(
"mixture").subDict(
"transport")
192 *2.0*gamma_/
Pr.
value()/(gamma_ + 1.0)
193 *(2.0 - accommodationCoeff_)/accommodationCoeff_
196 Field<scalar> aCoeff(prho.snGrad() - prho/C2);
197 Field<scalar> KEbyRho(0.5*
magSqr(pU));
199 valueFraction() = (1.0/(1.0 +
patch().deltaCoeffs()*C2));
203 mixedFvPatchScalarField::updateCoeffs();
217 os.
writeEntry(
"accommodationCoeff", accommodationCoeff_);
220 writeEntry(
"value",
os);
231 smoluchowskiJumpTFvPatchScalarField
Macros for easy insertion into run-time selection tables.
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.
static const word dictName
const Type & value() const
Return const reference to value.
virtual bool write()
Write the output fields.
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Smoluchowski temperature jump boundary condition.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
constexpr scalar piByTwo(0.5 *M_PI)
const std::string patch
OpenFOAM patch number as a std::string.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< vector, fvPatchField, volMesh > volVectorField
const dimensionSet dimless
Dimensionless.
List< label > labelList
A List of labels.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
fvPatchField< vector > fvPatchVectorField
fvPatchField< scalar > fvPatchScalarField
dimensionedScalar Pr("Pr", dimless, laminarTransport)