26 #include "JohnsonJacksonParticleThetaFvPatchScalarField.H" 28 #include "twoPhaseSystem.H" 37 JohnsonJacksonParticleThetaFvPatchScalarField
47 const DimensionedField<scalar, volMesh>& iF
50 mixedFvPatchScalarField(p, iF),
51 restitutionCoefficient_(
"restitutionCoefficient",
dimless, 0),
52 specularityCoefficient_(
"specularityCoefficient",
dimless, 0)
59 const JohnsonJacksonParticleThetaFvPatchScalarField& ptf,
61 const DimensionedField<scalar, volMesh>& iF,
62 const fvPatchFieldMapper& mapper
65 mixedFvPatchScalarField(ptf, p, iF, mapper),
66 restitutionCoefficient_(ptf.restitutionCoefficient_),
67 specularityCoefficient_(ptf.specularityCoefficient_)
75 const DimensionedField<scalar, volMesh>& iF,
76 const dictionary& dict
79 mixedFvPatchScalarField(p, iF),
80 restitutionCoefficient_
82 "restitutionCoefficient",
84 dict.
lookup(
"restitutionCoefficient")
86 specularityCoefficient_
88 "specularityCoefficient",
90 dict.
lookup(
"specularityCoefficient")
95 (restitutionCoefficient_.value() < 0)
96 || (restitutionCoefficient_.value() > 1)
100 <<
"The restitution coefficient has to be between 0 and 1" 106 (specularityCoefficient_.
value() < 0)
107 || (specularityCoefficient_.
value() > 1)
111 <<
"The specularity coefficient has to be between 0 and 1" 125 const JohnsonJacksonParticleThetaFvPatchScalarField& ptf
128 mixedFvPatchScalarField(ptf),
129 restitutionCoefficient_(ptf.restitutionCoefficient_),
130 specularityCoefficient_(ptf.specularityCoefficient_)
137 const JohnsonJacksonParticleThetaFvPatchScalarField& ptf,
138 const DimensionedField<scalar, volMesh>& iF
141 mixedFvPatchScalarField(ptf, iF),
142 restitutionCoefficient_(ptf.restitutionCoefficient_),
143 specularityCoefficient_(ptf.specularityCoefficient_)
151 const fvPatchFieldMapper& m
154 mixedFvPatchScalarField::autoMap(m);
164 mixedFvPatchScalarField::rmap(ptf, addr);
176 const twoPhaseSystem& fluid = db().lookupObject<twoPhaseSystem>
181 const phaseModel& phased
183 fluid.phase1().name() == internalField().group()
191 patch().lookupPatchField<volScalarField, scalar>
193 phased.volScalarField::name()
199 patch().lookupPatchField<volVectorField, vector>
207 patch().lookupPatchField<volScalarField, scalar>
215 patch().lookupPatchField<volScalarField, scalar>
229 .lookupObject<IOdictionary>
234 .subDict(
"kineticTheoryCoeffs")
239 if (restitutionCoefficient_.value() != 1.0)
243 *specularityCoefficient_.
value()
245 /(scalar(1) -
sqr(restitutionCoefficient_.value()));
247 this->refGrad() = 0.0;
254 *(scalar(1) -
sqr(restitutionCoefficient_.value()))
259 this->valueFraction() =
c/(
c + patch().deltaCoeffs());
266 this->refValue() = 0.0;
271 *specularityCoefficient_.
value()
278 this->valueFraction() = 0.0;
281 mixedFvPatchScalarField::updateCoeffs();
291 os.writeKeyword(
"restitutionCoefficient")
293 os.writeKeyword(
"specularityCoefficient")
295 writeEntry(
"value", os);
fvPatchField< vector > fvPatchVectorField
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
const Type & value() const
Return const reference to value.
virtual void write(Ostream &) const
Write.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Macros for easy insertion into run-time selection tables.
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
dimensionedScalar pos(const dimensionedScalar &ds)
JohnsonJacksonParticleThetaFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
stressControl lookup("compactNormalStress") >> compactNormalStress
fvPatchField< scalar > fvPatchScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
List< label > labelList
A List of labels.
errorManip< error > abort(error &err)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensionedScalar alphaMax(laminarTransport.lookup("alphaMax"))
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionedScalar c
Speed of light in a vacuum.
virtual void updateCoeffs()
Update the coefficients.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void write(Ostream &) const
Write.