37 JohnsonJacksonParticleThetaFvPatchScalarField
47 const DimensionedField<scalar, volMesh>& iF
50 mixedFvPatchScalarField(p, iF),
51 restitutionCoefficient_(
"restitutionCoefficient",
dimless, 0),
52 specularityCoefficient_(
"specularityCoefficient",
dimless, 0)
60 const DimensionedField<scalar, volMesh>& iF,
61 const dictionary& dict
64 mixedFvPatchScalarField(p, iF),
65 restitutionCoefficient_
67 "restitutionCoefficient",
69 dict.
lookup(
"restitutionCoefficient")
71 specularityCoefficient_
73 "specularityCoefficient",
75 dict.
lookup(
"specularityCoefficient")
80 (restitutionCoefficient_.value() < 0)
81 || (restitutionCoefficient_.value() > 1)
85 <<
"The restitution coefficient has to be between 0 and 1" 91 (specularityCoefficient_.
value() < 0)
92 || (specularityCoefficient_.
value() > 1)
96 <<
"The specularity coefficient has to be between 0 and 1" 110 const JohnsonJacksonParticleThetaFvPatchScalarField& ptf,
112 const DimensionedField<scalar, volMesh>& iF,
113 const fvPatchFieldMapper& mapper
116 mixedFvPatchScalarField(ptf, p, iF, mapper),
117 restitutionCoefficient_(ptf.restitutionCoefficient_),
118 specularityCoefficient_(ptf.specularityCoefficient_)
126 const JohnsonJacksonParticleThetaFvPatchScalarField& ptf,
127 const DimensionedField<scalar, volMesh>& iF
130 mixedFvPatchScalarField(ptf, iF),
131 restitutionCoefficient_(ptf.restitutionCoefficient_),
132 specularityCoefficient_(ptf.specularityCoefficient_)
146 const phaseSystem& fluid =
149 const phaseModel& phase
151 fluid.phases()[internalField().group()]
157 patch().lookupPatchField<volScalarField, scalar>
159 phase.volScalarField::name()
165 patch().lookupPatchField<volVectorField, vector>
173 patch().lookupPatchField<volScalarField, scalar>
181 patch().lookupPatchField<volScalarField, scalar>
190 if (restitutionCoefficient_.value() != 1.0)
194 *specularityCoefficient_.
value()
196 /(scalar(1) -
sqr(restitutionCoefficient_.value()));
198 this->refGrad() = 0.0;
205 *(scalar(1) -
sqr(restitutionCoefficient_.value()))
207 /
max(4*
kappa*phase.alphaMax(), small)
210 this->valueFraction() =
c/(
c + patch().deltaCoeffs());
217 this->refValue() = 0.0;
222 *specularityCoefficient_.
value()
227 /
max(6*
kappa*phase.alphaMax(), small);
229 this->valueFraction() = 0;
232 mixedFvPatchScalarField::updateCoeffs();
242 writeEntry(os,
"restitutionCoefficient", restitutionCoefficient_);
243 writeEntry(os,
"specularityCoefficient", specularityCoefficient_);
static const word propertiesName
Default name of the phase properties dictionary.
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
fvPatchField< vector > fvPatchVectorField
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimless
virtual void write(Ostream &) const
Write.
const dimensionedScalar c
Speed of light in a vacuum.
Macros for easy insertion into run-time selection tables.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
JohnsonJacksonParticleThetaFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
stressControl lookup("compactNormalStress") >> compactNormalStress
virtual void write(Ostream &) const
Write.
fvPatchField< scalar > fvPatchScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
const Type & value() const
Return const reference to value.
errorManip< error > abort(error &err)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensionedScalar pos0(const dimensionedScalar &ds)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
virtual void updateCoeffs()
Update the coefficients.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)