31 template<
class Type,
class Limiter>
42 template<
class Type,
class Limiter>
61 template<
class Type,
class Limiter>
68 const VolField<Type>& vsf,
72 const fvMesh& mesh = vsf.mesh();
75 tGrad = basicGradScheme_().calcGrad(vsf,
name);
85 const labelUList& neighbour = mesh.neighbour();
90 Field<Type> maxVsf(vsf.primitiveField());
91 Field<Type> minVsf(vsf.primitiveField());
95 label own = owner[facei];
96 label nei = neighbour[facei];
98 const Type& vsfOwn = vsf[own];
99 const Type& vsfNei = vsf[nei];
101 maxVsf[own] =
max(maxVsf[own], vsfNei);
102 minVsf[own] =
min(minVsf[own], vsfNei);
104 maxVsf[nei] =
max(maxVsf[nei], vsfOwn);
105 minVsf[nei] =
min(minVsf[nei], vsfOwn);
109 const typename VolField<Type>::Boundary& bsf =
114 const fvPatchField<Type>& psf = bsf[
patchi];
119 const Field<Type> psfNei(psf.patchNeighbourField());
123 label own = pOwner[pFacei];
124 const Type& vsfNei = psfNei[pFacei];
126 maxVsf[own] =
max(maxVsf[own], vsfNei);
127 minVsf[own] =
min(minVsf[own], vsfNei);
134 label own = pOwner[pFacei];
135 const Type& vsfNei = psf[pFacei];
137 maxVsf[own] =
max(maxVsf[own], vsfNei);
138 minVsf[own] =
min(minVsf[own], vsfNei);
148 const Field<Type> maxMinVsf((1.0/k_ - 1.0)*(maxVsf - minVsf));
156 Field<Type>
limiter(vsf.primitiveField().size(), pTraits<Type>::one);
160 label own = owner[facei];
161 label nei = neighbour[facei];
169 (Cf[facei] -
C[own]) & g[own]
178 (Cf[facei] -
C[nei]) & g[nei]
189 label own = pOwner[pFacei];
196 ((pCf[pFacei] -
C[own]) & g[own])
203 Info<<
"gradient limiter for: " << vsf.name()
210 g.correctBoundaryConditions();
static const Foam::dimensionedScalar C("C", Foam::dimTemperature, 234.5)
#define forAll(list, i)
Loop across all elements in list.
Generic GeometricField class.
cellLimitedGrad gradient scheme applied to a runTime selected base gradient scheme.
virtual tmp< VolField< typename outerProduct< vector, Type >::type > > calcGrad(const VolField< Type > &vsf, const word &name) const
Return the gradient of the given field to the gradScheme::grad.
A class for managing temporary objects.
volVectorField vectorField(fieldObject, mesh)
U correctBoundaryConditions()
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)
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.
Tensor< scalar > tensor
Tensor of scalars.
Ostream & endl(Ostream &os)
Add newline and flush stream.
word name(const bool)
Return a word representation of a bool.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
Type gAverage(const FieldField< Field, Type > &f)
Type gMin(const FieldField< Field, Type > &f)
SurfaceField< vector > surfaceVectorField
UList< label > labelUList
Type gMax(const FieldField< Field, Type > &f)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.