41 const DimensionedField<vector, volMesh>& iF
44 partialSlipFvPatchVectorField(
p, iF),
47 psiName_(
"thermo:psi"),
50 accommodationCoeff_(1.0),
51 Uwall_(
p.size(), Zero),
59 const maxwellSlipUFvPatchVectorField& mspvf,
61 const DimensionedField<vector, volMesh>& iF,
62 const fvPatchFieldMapper& mapper
65 partialSlipFvPatchVectorField(mspvf,
p, iF, mapper),
67 rhoName_(mspvf.rhoName_),
68 psiName_(mspvf.psiName_),
69 muName_(mspvf.muName_),
70 tauMCName_(mspvf.tauMCName_),
71 accommodationCoeff_(mspvf.accommodationCoeff_),
73 thermalCreep_(mspvf.thermalCreep_),
74 curvature_(mspvf.curvature_)
81 const DimensionedField<vector, volMesh>& iF,
82 const dictionary&
dict
85 partialSlipFvPatchVectorField(
p, iF),
86 TName_(
dict.getOrDefault<word>(
"T",
"T")),
87 rhoName_(
dict.getOrDefault<word>(
"rho",
"rho")),
88 psiName_(
dict.getOrDefault<word>(
"psi",
"thermo:psi")),
89 muName_(
dict.getOrDefault<word>(
"mu",
"thermo:mu")),
90 tauMCName_(
dict.getOrDefault<word>(
"tauMC",
"tauMC")),
91 accommodationCoeff_(
dict.
get<scalar>(
"accommodationCoeff")),
92 Uwall_(
"Uwall",
dict,
p.size()),
93 thermalCreep_(
dict.getOrDefault(
"thermalCreep", true)),
94 curvature_(
dict.getOrDefault(
"curvature", true))
98 mag(accommodationCoeff_) < SMALL
99 ||
mag(accommodationCoeff_) > 2.0
103 <<
"unphysical accommodationCoeff_ specified"
104 <<
"(0 < accommodationCoeff_ <= 1)" <<
endl
105 <<
exit(FatalIOError);
108 if (
dict.found(
"value"))
110 fvPatchField<vector>::operator=
115 if (
dict.found(
"refValue") &&
dict.found(
"valueFraction"))
118 this->valueFraction() =
123 this->refValue() = *
this;
124 this->valueFraction() = scalar(1);
132 const maxwellSlipUFvPatchVectorField& mspvf,
133 const DimensionedField<vector, volMesh>& iF
136 partialSlipFvPatchVectorField(mspvf, iF),
137 TName_(mspvf.TName_),
138 rhoName_(mspvf.rhoName_),
139 psiName_(mspvf.psiName_),
140 muName_(mspvf.muName_),
141 tauMCName_(mspvf.tauMCName_),
142 accommodationCoeff_(mspvf.accommodationCoeff_),
143 Uwall_(mspvf.Uwall_),
144 thermalCreep_(mspvf.thermalCreep_),
145 curvature_(mspvf.curvature_)
162 const fvPatchField<scalar>& ppsi =
168 * (2.0 - accommodationCoeff_)/accommodationCoeff_
171 Field<scalar> pnu(pmu/prho);
172 valueFraction() = (1.0/(1.0 +
patch().deltaCoeffs()*C1*pnu));
180 label patchi = this->
patch().index();
182 Field<vector> gradpT(
fvc::grad(vsfT)().boundaryField()[patchi]);
185 refValue() -= 3.0*pnu/(4.0*pT)*
transform(
I -
n*
n, gradpT);
197 partialSlipFvPatchVectorField::updateCoeffs();
210 os.
writeEntry(
"accommodationCoeff", accommodationCoeff_);
218 writeEntry(
"value",
os);
229 maxwellSlipUFvPatchVectorField
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.
virtual bool write()
Write the output fields.
Maxwell slip boundary condition including thermal creep and surface curvature terms that can be optio...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Tensor of scalars, i.e. Tensor<scalar>.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
Calculate the gradient of the given field.
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.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
static const Identity< scalar > I
fvPatchField< tensor > fvPatchTensorField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
GeometricField< tensor, fvPatchField, volMesh > volTensorField
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Field< vector > vectorField
Specialisation of Field<T> for vector.
errorManipArg< error, int > exit(error &err, const int errNo=1)
fvPatchField< vector > fvPatchVectorField
fvPatchField< scalar > fvPatchScalarField