37namespace kineticTheoryModels
39namespace frictionalStressModels
45 frictionalStressModel,
46 JohnsonJacksonSchaeffer,
56Foam::kineticTheoryModels::frictionalStressModels::
57JohnsonJacksonSchaeffer::JohnsonJacksonSchaeffer
59 const dictionary&
dict
62 frictionalStressModel(
dict),
63 coeffDict_(
dict.optionalSubDict(typeName +
"Coeffs")),
64 Fr_(
"Fr", dimensionSet(1, -1, -2, 0, 0), coeffDict_),
65 eta_(
"eta",
dimless, coeffDict_),
67 phi_(
"phi",
dimless, coeffDict_),
68 alphaDeltaMin_(
"alphaDeltaMin",
dimless, coeffDict_)
87 const phaseModel& phase,
88 const dimensionedScalar& alphaMinFriction,
95 Fr_*
pow(
max(
alpha - alphaMinFriction, scalar(0)), eta_)
104 const phaseModel& phase,
105 const dimensionedScalar& alphaMinFriction,
113 eta_*
pow(
max(
alpha - alphaMinFriction, scalar(0)), eta_ - 1.0)
115 + p_*
pow(
max(
alpha - alphaMinFriction, scalar(0)), eta_)
124 const phaseModel& phase,
125 const dimensionedScalar& alphaMinFriction,
127 const volScalarField& pf,
128 const volSymmTensorField&
D
133 tmp<volScalarField> tnu
139 "JohnsonJacksonSchaeffer:nu",
140 phase.mesh().time().timeName(),
155 if (
alpha[celli] > alphaMinFriction.value())
158 0.5*pf[celli]*
sin(phi_.value())
169 volScalarField::Boundary& nufBf = nuf.boundaryFieldRef();
177 pf.boundaryField()[patchi]*
sin(phi_.value())
179 mag(
U.boundaryField()[patchi].snGrad())
187 nuf.correctBoundaryConditions();
196 coeffDict_ <<= dict_.optionalSubDict(typeName +
"Coeffs");
198 Fr_.read(coeffDict_);
199 eta_.read(coeffDict_);
202 phi_.read(coeffDict_);
205 alphaDeltaMin_.read(coeffDict_);
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
virtual tmp< volScalarField > frictionalPressure(const phaseModel &phase, const dimensionedScalar &alphaMinFriction, const dimensionedScalar &alphaMax) const
virtual tmp< volScalarField > frictionalPressurePrime(const phaseModel &phase, const dimensionedScalar &alphaMinFriction, const dimensionedScalar &alphaMax) const
virtual ~JohnsonJacksonSchaeffer()
Destructor.
A class for managing temporary objects.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
const polyBoundaryMesh & patches
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< vector, fvPatchField, volMesh > volVectorField
const dimensionSet dimless
Dimensionless.
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
PtrList< fvPatch > fvPatchList
Store lists of fvPatch as a PtrList.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Cmpt invariantII(const SymmTensor< Cmpt > &st)
Return the 2nd invariant of a SymmTensor.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar & D
dimensionedScalar alphaMax("alphaMax", dimless/dimTime, laminarTransport)
#define forAll(list, i)
Loop across all elements in list.
Unit conversion functions.