34 template<
class Type,
class Limiter>
45 template<
class Type,
class Limiter>
64 template<
class Type,
class Limiter>
76 const GeometricField<Type, fvPatchField, volMesh>& vsf,
80 const fvMesh&
mesh = vsf.mesh();
86 > tGrad = basicGradScheme_().calcGrad(vsf,
name);
106 Field<Type> maxVsf(vsf.primitiveField());
107 Field<Type> minVsf(vsf.primitiveField());
111 const label own = owner[facei];
112 const label nei = neighbour[facei];
114 const Type& vsfOwn = vsf[own];
115 const Type& vsfNei = vsf[nei];
117 maxVsf[own] =
max(maxVsf[own], vsfNei);
118 minVsf[own] =
min(minVsf[own], vsfNei);
120 maxVsf[nei] =
max(maxVsf[nei], vsfOwn);
121 minVsf[nei] =
min(minVsf[nei], vsfOwn);
125 const auto& bsf = vsf.boundaryField();
129 const fvPatchField<Type>& psf = bsf[patchi];
134 const Field<Type> psfNei(psf.patchNeighbourField());
138 const label own = pOwner[pFacei];
139 const Type& vsfNei = psfNei[pFacei];
141 maxVsf[own] =
max(maxVsf[own], vsfNei);
142 minVsf[own] =
min(minVsf[own], vsfNei);
149 const label own = pOwner[pFacei];
150 const Type& vsfNei = psf[pFacei];
152 maxVsf[own] =
max(maxVsf[own], vsfNei);
153 minVsf[own] =
min(minVsf[own], vsfNei);
163 const Field<Type> maxMinVsf((1.0/k_ - 1.0)*(maxVsf - minVsf));
171 Field<Type>
limiter(vsf.primitiveField().size(), pTraits<Type>::one);
175 const label own = owner[facei];
176 const label nei = neighbour[facei];
184 (Cf[facei] -
C[own]) &
g[own]
193 (Cf[facei] -
C[nei]) &
g[nei]
200 const vectorField& pCf = Cf.boundaryField()[patchi];
204 const label own = pOwner[pFacei];
211 ((pCf[pFacei] -
C[own]) &
g[own])
218 Info<<
"gradient limiter for: " << vsf.name()
225 g.correctBoundaryConditions();
virtual tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > calcGrad(const GeometricField< Type, fvPatchField, volMesh > &vsf, const word &name) const
Return the gradient of the given field to the gradScheme::grad.
Type gMin(const FieldField< Field, Type > &f)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
cellMask correctBoundaryConditions()
Ostream & endl(Ostream &os)
Add newline and flush stream.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
typeOfRank< typename pTraits< arg1 >::cmptType, direction(pTraits< arg1 >::rank)+direction(pTraits< arg2 >::rank) >::type type
Generic GeometricField class.
tmp< areaScalarField > limiter(const areaScalarField &phi)
GeometricField< vector, fvPatchField, volMesh > volVectorField
UList< label > labelUList
A UList of labels.
#define forAll(list, i)
Loop across all elements in list.
Mesh data needed to do the Finite Volume discretisation.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
cellLimitedGrad gradient scheme applied to a runTime selected base gradient scheme.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
const uniformDimensionedVectorField & g
int debug
Static debugging option.
Type gMax(const FieldField< Field, Type > &f)
Type gAverage(const FieldField< Field, Type > &f)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< vector > vectorField
Specialisation of Field<T> for vector.
A class for managing temporary objects.