33 template<
class Type,
class Limiter,
template<
class>
class LimitFunc>
36 const VolField<Type>& phi,
40 const fvMesh& mesh = this->mesh();
42 tmp<VolField<typename Limiter::phiType>> tlPhi = LimitFunc<Type>()(phi);
43 const VolField<typename Limiter::phiType>& lPhi = tlPhi();
45 tmp<VolField<typename Limiter::gradPhiType>> tgradc(
fvc::grad(lPhi));
46 const VolField<typename Limiter::gradPhiType>& gradc = tgradc();
51 const labelUList& neighbour = mesh.neighbour();
55 scalarField& pLim = limiterField.primitiveFieldRef();
59 label own = owner[face];
60 label nei = neighbour[face];
65 this->faceFlux_[face],
74 const typename VolField<Type>::Boundary&
75 bPhi = phi.boundaryField();
77 surfaceScalarField::Boundary& bLim =
78 limiterField.boundaryFieldRef();
84 if (bPhi[
patchi].coupled())
88 this->faceFlux_.boundaryField()[
patchi];
90 const Field<typename Limiter::phiType> plPhiP
92 lPhi.boundaryField()[
patchi].patchInternalField()
94 const Field<typename Limiter::phiType> plPhiN
96 lPhi.boundaryField()[
patchi].patchNeighbourField()
98 const Field<typename Limiter::gradPhiType> pGradcP
100 gradc.boundaryField()[
patchi].patchInternalField()
102 const Field<typename Limiter::gradPhiType> pGradcN
104 gradc.boundaryField()[
patchi].patchNeighbourField()
134 template<
class Type,
class Limiter,
template<
class>
class LimitFunc>
141 const fvMesh& mesh = this->mesh();
143 const word limiterFieldName(
type() +
"Limiter(" + phi.
name() +
')');
145 if (this->mesh().
solution().cache(
"limiter"))
166 mesh.objectRegistry::store(limiterField);
175 calcLimiter(phi, limiterField);
191 calcLimiter(phi, tlimiterField.
ref());
193 return tlimiterField;
static const Foam::dimensionedScalar C("C", Foam::dimTemperature, 234.5)
#define forAll(list, i)
Loop across all elements in list.
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
const word & name() const
Return name.
Class to create NVD/TVD limited weighting-factors.
virtual tmp< surfaceScalarField > limiter(const VolField< Type > &) const
Return the interpolation weighting factors.
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.
Type & lookupObjectRef(const word &name) const
Lookup and return the object reference of the given Type.
bool foundObject(const word &name) const
Is the named Type in registry.
Selector class for relaxation factors, solver type and solution.
A class for managing temporary objects.
T & ref() const
Return non-const reference or generate a fatal error.
A class for handling words, derived from string.
Calculate the gradient of the given field.
volScalarField scalarField(fieldObject, mesh)
volVectorField vectorField(fieldObject, mesh)
void limiter(surfaceScalarField &lambda, const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phiBD, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su, const PsiMaxType &psiMax, const PsiMinType &psiMin)
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const dimensionSet dimless
SurfaceField< scalar > surfaceScalarField
UList< label > labelUList
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.