30 #include "twoPhaseSystem.H" 39 JohnsonJacksonParticleThetaFvPatchScalarField
52 mixedFvPatchScalarField(
p, iF),
53 restitutionCoefficient_(
"restitutionCoefficient",
dimless,
Zero),
54 specularityCoefficient_(
"specularityCoefficient",
dimless,
Zero)
67 mixedFvPatchScalarField(ptf,
p, iF, mapper),
68 restitutionCoefficient_(ptf.restitutionCoefficient_),
69 specularityCoefficient_(ptf.specularityCoefficient_)
82 mixedFvPatchScalarField(
p, iF),
83 restitutionCoefficient_(
"restitutionCoefficient",
dimless,
dict),
84 specularityCoefficient_(
"specularityCoefficient",
dimless,
dict)
88 (restitutionCoefficient_.
value() < 0)
89 || (restitutionCoefficient_.
value() > 1)
93 <<
"The restitution coefficient has to be between 0 and 1" 99 (specularityCoefficient_.
value() < 0)
100 || (specularityCoefficient_.
value() > 1)
104 <<
"The specularity coefficient has to be between 0 and 1" 118 mixedFvPatchScalarField(ptf),
119 restitutionCoefficient_(ptf.restitutionCoefficient_),
120 specularityCoefficient_(ptf.specularityCoefficient_)
131 mixedFvPatchScalarField(ptf, iF),
132 restitutionCoefficient_(ptf.restitutionCoefficient_),
133 specularityCoefficient_(ptf.specularityCoefficient_)
144 mixedFvPatchScalarField::autoMap(m);
154 mixedFvPatchScalarField::rmap(ptf, addr);
166 const twoPhaseSystem&
fluid = db().lookupObject<twoPhaseSystem>
171 const phaseModel& phased
181 patch().lookupPatchField<volScalarField>
183 phased.volScalarField::name()
189 patch().lookupPatchField<volVectorField>
197 patch().lookupPatchField<volScalarField>
205 patch().lookupPatchField<volScalarField>
219 .lookupObject<IOdictionary>
224 .subDict(
"kineticTheoryCoeffs")
228 if (restitutionCoefficient_.value() != 1.0)
232 *specularityCoefficient_.
value()
234 /(scalar(1) -
sqr(restitutionCoefficient_.value()));
236 this->refGrad() = 0.0;
243 *(scalar(1) -
sqr(restitutionCoefficient_.value()))
248 this->valueFraction() =
c/(
c +
patch().deltaCoeffs());
255 this->refValue() = 0.0;
260 *specularityCoefficient_.value()
267 this->valueFraction() = 0;
270 mixedFvPatchScalarField::updateCoeffs();
280 os.
writeEntry(
"restitutionCoefficient", restitutionCoefficient_);
281 os.
writeEntry(
"specularityCoefficient", specularityCoefficient_);
const phaseModel & phase1() const
Return phase model 1.
const Type & value() const noexcept
Return const reference to value.
fvPatchField< vector > fvPatchVectorField
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
void writeValueEntry(Ostream &os) const
Write *this field as a "value" entry.
const dimensionSet dimless
Dimensionless.
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.
const phaseModel & phase2() const
Return phase model 2.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
dimensionedScalar alphaMax("alphaMax", dimless/dimTime, laminarTransport)
JohnsonJacksonParticleThetaFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
fvPatchField< scalar > fvPatchScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
A FieldMapper for finite-volume patch fields.
constexpr scalar pi(M_PI)
errorManip< error > abort(error &err)
dimensionedScalar pos0(const dimensionedScalar &ds)
OBJstream os(runTime.globalPath()/outputName)
virtual void write(Ostream &) const
Write.
virtual void updateCoeffs()
Update the coefficients.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionedScalar c
Speed of light in a vacuum.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const std::string patch
OpenFOAM patch number as a std::string.
const word & name() const
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Robin condition for the particulate granular temperature.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)