47 const volScalarField& vsf,
51 const fvMesh&
mesh = vsf.mesh();
53 tmp<volVectorField> tGrad = basicGradScheme_().calcGrad(vsf, name);
60 volVectorField&
g = tGrad.ref();
62 const labelUList& owner =
mesh.owner();
63 const labelUList& neighbour =
mesh.neighbour();
65 const volVectorField&
C =
mesh.C();
66 const surfaceVectorField& Cf =
mesh.Cf();
69 scalarField limiter(vsf.primitiveField().size(), 1.0);
71 const scalar rk = (1.0/k_ - 1.0);
75 const label own = owner[facei];
76 const label nei = neighbour[facei];
78 const scalar vsfOwn = vsf[own];
79 const scalar vsfNei = vsf[nei];
81 scalar maxFace =
max(vsfOwn, vsfNei);
82 scalar minFace = min(vsfOwn, vsfNei);
83 const scalar maxMinFace = rk*(maxFace - minFace);
84 maxFace += maxMinFace;
85 minFace -= maxMinFace;
91 maxFace - vsfOwn, minFace - vsfOwn,
92 (Cf[facei] -
C[own]) &
g[own]
99 maxFace - vsfNei, minFace - vsfNei,
100 (Cf[facei] -
C[nei]) &
g[nei]
104 const volScalarField::Boundary& bsf = vsf.boundaryField();
108 const fvPatchScalarField& psf = bsf[patchi];
110 const labelUList& pOwner =
mesh.boundary()[patchi].faceCells();
111 const vectorField& pCf = Cf.boundaryField()[patchi];
115 const scalarField psfNei(psf.patchNeighbourField());
119 const label own = pOwner[pFacei];
121 const scalar vsfOwn = vsf[own];
122 const scalar vsfNei = psfNei[pFacei];
124 scalar maxFace =
max(vsfOwn, vsfNei);
125 scalar minFace = min(vsfOwn, vsfNei);
126 const scalar maxMinFace = rk*(maxFace - minFace);
127 maxFace += maxMinFace;
128 minFace -= maxMinFace;
133 maxFace - vsfOwn, minFace - vsfOwn,
134 (pCf[pFacei] -
C[own]) &
g[own]
138 else if (psf.fixesValue())
142 const label own = pOwner[pFacei];
144 const scalar vsfOwn = vsf[own];
145 const scalar vsfNei = psf[pFacei];
147 scalar maxFace =
max(vsfOwn, vsfNei);
148 scalar minFace = min(vsfOwn, vsfNei);
149 const scalar maxMinFace = rk*(maxFace - minFace);
150 maxFace += maxMinFace;
151 minFace -= maxMinFace;
156 maxFace - vsfOwn, minFace - vsfOwn,
157 (pCf[pFacei] -
C[own]) &
g[own]
165 Info<<
"gradient limiter for: " << vsf.name()
166 <<
" max = " << gMax(limiter)
167 <<
" min = " << gMin(limiter)
168 <<
" average: " << gAverage(limiter) << endl;
171 g.primitiveFieldRef() *= limiter;
172 g.correctBoundaryConditions();
173 gaussGrad<scalar>::correctBoundaryConditions(vsf,
g);
207 const scalar rk = (1.0/k_ - 1.0);
211 const label own = owner[facei];
212 const label nei = neighbour[facei];
214 const vector& vvfOwn = vvf[own];
215 const vector& vvfNei = vvf[nei];
218 vector gradf((Cf[facei] -
C[own]) &
g[own]);
220 scalar vsfOwn = gradf & vvfOwn;
221 scalar vsfNei = gradf & vvfNei;
223 scalar maxFace =
max(vsfOwn, vsfNei);
224 scalar minFace =
min(vsfOwn, vsfNei);
225 const scalar maxMinFace = rk*(maxFace - minFace);
226 maxFace += maxMinFace;
227 minFace -= maxMinFace;
232 maxFace - vsfOwn, minFace - vsfOwn,
238 gradf = (Cf[facei] -
C[nei]) &
g[nei];
240 vsfOwn = gradf & vvfOwn;
241 vsfNei = gradf & vvfNei;
243 maxFace =
max(vsfOwn, vsfNei);
244 minFace =
min(vsfOwn, vsfNei);
249 maxFace - vsfNei, minFace - vsfNei,
270 const label own = pOwner[pFacei];
272 const vector& vvfOwn = vvf[own];
273 const vector& vvfNei = psfNei[pFacei];
275 const vector gradf((pCf[pFacei] -
C[own]) &
g[own]);
277 const scalar vsfOwn = gradf & vvfOwn;
278 const scalar vsfNei = gradf & vvfNei;
280 scalar maxFace =
max(vsfOwn, vsfNei);
281 scalar minFace =
min(vsfOwn, vsfNei);
282 const scalar maxMinFace = rk*(maxFace - minFace);
283 maxFace += maxMinFace;
284 minFace -= maxMinFace;
289 maxFace - vsfOwn, minFace - vsfOwn,
298 const label own = pOwner[pFacei];
300 const vector& vvfOwn = vvf[own];
301 const vector& vvfNei = psf[pFacei];
303 const vector gradf((pCf[pFacei] -
C[own]) &
g[own]);
305 const scalar vsfOwn = gradf & vvfOwn;
306 const scalar vsfNei = gradf & vvfNei;
308 scalar maxFace =
max(vsfOwn, vsfNei);
309 scalar minFace =
min(vsfOwn, vsfNei);
310 const scalar maxMinFace = rk*(maxFace - minFace);
311 maxFace += maxMinFace;
312 minFace -= maxMinFace;
317 maxFace - vsfOwn, minFace - vsfOwn,
326 Info<<
"gradient limiter for: " << vvf.
name()
333 g.correctBoundaryConditions();
static const Foam::dimensionedScalar C("", Foam::dimTemperature, 234.5)
const uniformDimensionedVectorField & g
Graphite solid properties.
const Mesh & mesh() const
Return mesh.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const word & name() const noexcept
Return the object name.
void size(const label n)
Older name for setAddressableSize.
Mesh data needed to do the Finite Volume discretisation.
virtual bool fixesValue() const
Return true if this patch field fixes a value.
virtual tmp< Field< Type > > patchNeighbourField() const
Return patchField on the opposite patch of a coupled patch.
virtual bool coupled() const
Return true if this patch field is coupled.
faceLimitedGrad gradient scheme applied to a runTime selected base gradient scheme.
Basic second-order gradient scheme using face-interpolation and Gauss' theorem.
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define makeFvGradScheme(SS)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Type gAverage(const FieldField< Field, Type > &f)
Type gMin(const FieldField< Field, Type > &f)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Type gMax(const FieldField< Field, Type > &f)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
tmp< areaScalarField > limiter(const areaScalarField &phi)
#define forAll(list, i)
Loop across all elements in list.