41 template<
class BasePhaseModel>
45 word phiName(IOobject::groupName(
"phi", this->
name()));
47 typeIOobject<surfaceScalarField> phiHeader
50 U.mesh().time().name(),
55 if (phiHeader.headerOk())
57 Info<<
"Reading face flux field " << phiName <<
endl;
59 return tmp<surfaceScalarField>
66 U.mesh().time().name(),
77 Info<<
"Calculating face flux field " << phiName <<
endl;
81 U.boundaryField().size(),
82 calculatedFvPatchScalarField::typeName
87 if (!
U.boundaryField()[
patchi].assignable())
89 phiTypes[
patchi] = fixedValueFvPatchScalarField::typeName;
93 return tmp<surfaceScalarField>
100 U.mesh().time().name(),
115 template<
class BasePhaseModel>
119 const word& phaseName,
120 const bool referencePhase,
124 BasePhaseModel(fluid, phaseName, referencePhase, index),
130 fluid.mesh().time().
name(),
143 fluid.mesh().time().
name(),
156 fluid.mesh().time().
name(),
180 thermophysicalTransport_
186 >::
New(momentumTransport_, this->thermo_)
193 fluid.mesh().time().
name(),
203 if (fluid.
mesh().
dynamic() || this->fluid().MRF().size())
225 template<
class BasePhaseModel>
232 template<
class BasePhaseModel>
243 template<
class BasePhaseModel>
250 template<
class BasePhaseModel>
253 BasePhaseModel::correctKinematics();
269 K_.ref() = 0.5*
magSqr(this->
U());
274 template<
class BasePhaseModel>
277 BasePhaseModel::predictMomentumTransport();
278 momentumTransport_->predict();
282 template<
class BasePhaseModel>
285 BasePhaseModel::predictThermophysicalTransport();
286 thermophysicalTransport_->predict();
290 template<
class BasePhaseModel>
293 BasePhaseModel::correctMomentumTransport();
294 momentumTransport_->correct();
298 template<
class BasePhaseModel>
301 BasePhaseModel::correctThermophysicalTransport();
302 thermophysicalTransport_->correct();
306 template<
class BasePhaseModel>
309 const fvMesh& mesh = this->fluid().mesh();
325 template<
class BasePhaseModel>
332 template<
class BasePhaseModel>
343 +
fvm::SuSp(-this->continuityError(), U_)
345 + momentumTransport_->divDevTau(U_)
350 template<
class BasePhaseModel>
364 + momentumTransport_->divDevTau(U_)
369 template<
class BasePhaseModel>
377 template<
class BasePhaseModel>
385 template<
class BasePhaseModel>
393 template<
class BasePhaseModel>
401 template<
class BasePhaseModel>
409 template<
class BasePhaseModel>
417 template<
class BasePhaseModel>
425 template<
class BasePhaseModel>
436 <<
"Uf has not been allocated."
444 template<
class BasePhaseModel>
455 <<
"Uf has not been allocated."
463 template<
class BasePhaseModel>
471 template<
class BasePhaseModel>
479 template<
class BasePhaseModel>
487 template<
class BasePhaseModel>
495 template<
class BasePhaseModel>
503 template<
class BasePhaseModel>
511 template<
class BasePhaseModel>
531 template<
class BasePhaseModel>
541 byDt(phi_ - phi_.oldTime())
549 template<
class BasePhaseModel>
553 return continuityError_;
557 template<
class BasePhaseModel>
575 template<
class BasePhaseModel>
583 template<
class BasePhaseModel>
599 template<
class BasePhaseModel>
603 return momentumTransport_->k();
607 template<
class BasePhaseModel>
611 return momentumTransport_->pPrime();
615 template<
class BasePhaseModel>
619 return thermophysicalTransport_->kappaEff(
patchi);
623 template<
class BasePhaseModel>
627 return thermophysicalTransport_->divq(
he);
631 template<
class BasePhaseModel>
635 return thermophysicalTransport_->divj(Yi);
#define forAll(list, i)
Loop across all elements in list.
Generic GeometricField class.
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
writeOption writeOpt() const
static word groupName(Name name, const word &group)
virtual void correctUf()
Correct the face velocity for moving meshes.
virtual void correctKinematics()
Correct the kinematics.
virtual tmp< fvVectorMatrix > UEqn()
Return the momentum equation.
virtual tmp< surfaceScalarField > DUDtf() const
Return the substantive acceleration on the faces.
MovingPhaseModelTransportThermoModel< typename BasePhaseModel::thermoModel >::type transportThermoModel
Thermo type for the thermophysical transport model.
volVectorField U_
Velocity field.
virtual void correct()
Correct the phase properties other than the thermo.
virtual tmp< volScalarField > k() const
Return the turbulent kinetic energy.
virtual void predictThermophysicalTransport()
Predict the energy transport e.g. alphat.
virtual ~MovingPhaseModel()
Destructor.
virtual tmp< volVectorField > U() const
Return the velocity.
virtual void correctMomentumTransport()
Correct the momentumTransport.
virtual tmp< fvScalarMatrix > divj(volScalarField &Yi) const
Return the source term for the given specie mass-fraction.
virtual surfaceScalarField & alphaRhoPhiRef()
Access the mass flux of the phase.
virtual const autoPtr< volScalarField > & divU() const
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
virtual void correctThermophysicalTransport()
Correct the energy transport e.g. alphat.
virtual surfaceScalarField & phiRef()
Access the volumetric flux.
virtual tmp< surfaceScalarField > alphaPhi() const
Return the volumetric flux of the phase.
virtual tmp< volScalarField > K() const
Return the phase kinetic energy.
virtual tmp< volScalarField > pPrime() const
Return the phase-pressure'.
surfaceScalarField phi_
Flux.
MovingPhaseModel(const phaseSystem &fluid, const word &phaseName, const bool referencePhase, const label index)
virtual const autoPtr< surfaceVectorField > & Uf() const
Return the face velocity.
virtual surfaceVectorField & UfRef()
Access the face velocity.
virtual surfaceScalarField & alphaPhiRef()
Access the volumetric flux of the phase.
virtual tmp< surfaceScalarField > alphaRhoPhi() const
Return the mass flux of the phase.
virtual tmp< surfaceScalarField > phi() const
Return the volumetric flux.
virtual tmp< fvVectorMatrix > UfEqn()
Return the momentum equation for the face-based algorithm.
virtual volVectorField & URef()
Access the velocity.
virtual tmp< volVectorField > DUDt() const
Return the substantive acceleration.
virtual bool stationary() const
Return whether the phase is stationary.
virtual void predictMomentumTransport()
Predict the momentumTransport.
virtual tmp< volScalarField > continuityError() const
Return the continuity error.
virtual tmp< scalarField > kappaEff(const label patchi) const
Return the effective thermal conductivity on a patch.
virtual void correctContinuityError(const volScalarField &source)
Correct the continuity error.
autoPtr< surfaceVectorField > Uf_
Face velocity field.
virtual tmp< fvScalarMatrix > divq(volScalarField &he) const
Return the source term for the energy equation.
Templated base class for multiphase thermophysical transport models.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Dimension set for the base types.
const word & name() const
Return const reference to name.
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return the top-level database.
const surfaceVectorField & Sf() const
Return cell face area vectors.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
bool dynamic() const
Is this mesh dynamic?
Abstract base class for turbulence models (RAS, LES and laminar).
Class to represent a system of phases and model interfacial transfers between them.
const fvMesh & mesh() const
Return the mesh.
A class for managing temporary objects.
T * ptr() const
Return tmp pointer for reuse.
A class for handling words, derived from string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Calculate the first temporal derivative.
Calculate the divergence of the given field.
Calculate the face-flux of the given field.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for implicit and explicit sources.
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
tmp< SurfaceField< typename innerProduct< vector, Type >::type > > flux(const VolField< Type > &vf)
Return the face-flux field obtained from the given volVectorField.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< Type > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
tmp< VolField< Type > > DDt(const surfaceScalarField &phi, const VolField< Type > &psi)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
errorManipArg< error, int > exit(error &err, const int errNo=1)
List< word > wordList
A List of words.
VolField< vector > volVectorField
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
SurfaceField< scalar > surfaceScalarField
tmp< volScalarField > byDt(const volScalarField &vf)
const dimensionSet dimTime
const dimensionSet dimDensity
VolField< scalar > volScalarField
SurfaceField< vector > surfaceVectorField
word name(const complex &)
Return a string representation of a complex.
dimensioned< scalar > magSqr(const dimensioned< Type > &)