34 Foam::RASModels::kineticTheoryModel::continuousPhase()
const 36 const phaseSystem& fluid = phase_.
fluid();
40 if (fluid.movingPhases().size() != 2)
43 <<
"Continuous phase name must be specified " 44 <<
"when there are more than two moving phases." 48 forAll(fluid.movingPhases(), movingPhasei)
50 const phaseModel& otherPhase = fluid.movingPhases()[movingPhasei];
52 if (&otherPhase != &phase_)
59 return fluid.phases()[continuousPhaseName_];
96 kineticTheoryModels::viscosityModel::
New 103 kineticTheoryModels::conductivityModel::
New 110 kineticTheoryModels::radialModel::
New 115 granularPressureModel_
117 kineticTheoryModels::granularPressureModel::
New 122 frictionalStressModel_
124 kineticTheoryModels::frictionalStressModel::
New 149 dimensionSet(0, 2, -1, 0, 0),
222 if (type == typeName)
241 eddyViscosity<RASModel<phaseCompressibleMomentumTransportModel>>::
read()
244 coeffDict().lookup(
"equilibrium") >> equilibrium_;
245 e_.readIfPresent(coeffDict());
246 alphaMax_.readIfPresent(coeffDict());
247 alphaMinFriction_.readIfPresent(coeffDict());
249 viscosityModel_->read();
250 conductivityModel_->read();
251 radialModel_->read();
252 granularPressureModel_->read();
253 frictionalStressModel_->read();
283 return tmp<volSymmTensorField>
307 tmp<volScalarField> tpPrime
310 *granularPressureModel_->granularPressureCoeffPrime
313 radialModel_->g0(alpha_, alphaMinFriction_, alphaMax_),
314 radialModel_->g0prime(alpha_, alphaMinFriction_, alphaMax_),
318 + frictionalStressModel_->frictionalPressurePrime
326 volScalarField::Boundary& bpPrime =
327 tpPrime.
ref().boundaryFieldRef();
331 if (!bpPrime[
patchi].coupled())
351 return tmp<volSymmTensorField>
394 const phaseSystem& fluid = phase_.fluid();
395 const phaseModel& continuousPhase = this->continuousPhase();
405 tmp<volScalarField> tda(phase_.d());
408 tmp<volTensorField> tgradU(
fvc::grad(U_));
413 gs0_ = radialModel_->g0(alpha, alphaMinFriction_, alphaMax_);
418 nut_ = viscosityModel_->nu(alpha, Theta_, gs0_, rho, da, e_);
423 lambda_ = (4.0/3.0)*
sqr(alpha)*da*gs0_*(1 + e_)*ThetaSqrt/sqrtPi;
428 rho*(2*nut_*D + (lambda_ - (2.0/3.0)*nut_)*
tr(D)*I)
436 *
max(
sqr(alpha), residualAlpha_)
437 *rho*gs0_*(1.0/da)*ThetaSqrt/sqrtPi
443 fluid.foundSubModel<dragModel>(phase_, continuousPhase)
444 ? fluid.lookupSubModel<dragModel>(phase_, continuousPhase).K()
460 max(alpha, residualAlpha_)*rho
461 *sqrtPi*(ThetaSqrt + ThetaSmallSqrt)
468 granularPressureModel_->granularPressureCoeff
478 kappa_ = conductivityModel_->kappa(alpha, Theta_, gs0_, rho, da, e_);
496 -
fvm::SuSp((PsCoeff*I) && gradU, Theta_)
500 +
fvm::Sp(J2/(Theta_ + ThetaSmall), Theta_)
519 (sqrtPi/(3*(3.0 - e_)))
520 *(1 + 0.4*(1 + e_)*(3*e_ - 1)*alpha*gs0_)
521 +1.6*alpha*gs0_*(1 + e_)/sqrtPi
528 4*da*rho*(1 + e_)*alpha*gs0_/(3*sqrtPi) - 2*
K3/3.0
536 alpha/(alpha + residualAlpha_)
551 *(2*
K3*trD2 +
K2*tr2D)
557 /(2*
max(alpha, residualAlpha_)*
K4)
560 kappa_ = conductivityModel_->kappa(alpha, Theta_, gs0_, rho, da, e_);
568 nut_ = viscosityModel_->nu(alpha, Theta_, gs0_, rho, da, e_);
573 lambda_ = (4.0/3.0)*
sqr(alpha)*da*gs0_*(1 + e_)*ThetaSqrt/sqrtPi;
578 frictionalStressModel_->frictionalPressure
586 nuFric_ = frictionalStressModel_->nu
597 nuFric_ =
min(nuFric_, maxNut_ - nut_);
603 Info<< typeName <<
':' <<
nl 604 <<
" max(Theta) = " <<
max(Theta_).value() <<
nl 605 <<
" max(nut) = " <<
max(nut_).value() <<
endl;
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
virtual tmp< volSymmTensorField > devTau() const
Return the effective stress tensor.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
fvMatrix< scalar > fvScalarMatrix
#define forAll(list, i)
Loop across all elements in list.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
PhaseCompressibleMomentumTransportModel< phaseModel > phaseCompressibleMomentumTransportModel
Typedef for phaseCompressibleMomentumTransportModel.
static autoPtr< RASModel > New(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport)
Return a reference to the selected RAS model.
errorManipArg< error, int > exit(error &err, const int errNo=1)
virtual tmp< volSymmTensorField > sigma() const
Return the stress tensor [m^2/s^2].
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< surfaceScalarField > pPrimef() const
Return the face-phase-pressure'.
GeometricField< tensor, fvPatchField, volMesh > volTensorField
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
const Time & time() const
static const dictionary null
Null dictionary.
virtual tmp< fvVectorMatrix > divDevTau(volVectorField &U) const
Return the source term for the momentum equation.
IOobject(const word &name, const fileName &instance, const objectRegistry ®istry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< GeometricField< Type, fvPatchField, volMesh > > Sp(const volScalarField &sp, const GeometricField< Type, fvPatchField, volMesh > &vf)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
GeometricField< vector, fvPatchField, volMesh > volVectorField
const fvMesh & mesh() const
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dictionary coeffDict_
Model coefficients dictionary.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
RASModel< phaseCompressibleMomentumTransportModel > ::transportModel transportModel
static const SymmTensor I
static word groupName(Name name, const word &group)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
static const word null
An empty word.
virtual void correct()
Solve the kinetic theory equations and correct the viscosity.
void min(const dimensioned< Type > &)
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
const word & name() const
Name function is needed to disambiguate those inherited.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual ~kineticTheoryModel()
Destructor.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Internal & ref()
Return a reference to the dimensioned internal field.
kineticTheoryModel(const volScalarField &alpha, const volScalarField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const phaseModel &transport, const word &type=typeName)
Construct from components.
dimensionedSymmTensor dev2(const dimensionedSymmTensor &dt)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
static const dimensionSet dimK
Coefficient dimensions.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
virtual tmp< volScalarField > pPrime() const
Return the phase-pressure'.
const phaseSystem & fluid() const
Return the system to which this phase belongs.
eddyViscosity(const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport)
Construct from components.
virtual bool read()
Re-read model coefficients if they have changed.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
RASModel(const word &type, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport)
Construct from components.
const doubleScalar e
Elementary charge.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.