faceMDLimitedGrads.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "faceMDLimitedGrad.H"
27 #include "cellMDLimitedGrad.H"
28 #include "gaussGrad.H"
29 #include "fvMesh.H"
30 #include "volMesh.H"
31 #include "surfaceMesh.H"
32 #include "volFields.H"
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 makeFvGradScheme(faceMDLimitedGrad)
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 template<>
44 (
45  const volScalarField& vsf,
46  const word& name
47 ) const
48 {
49  const fvMesh& mesh = vsf.mesh();
50 
51  tmp<volVectorField> tGrad = basicGradScheme_().calcGrad(vsf, name);
52 
53  if (k_ < small)
54  {
55  return tGrad;
56  }
57 
58  volVectorField& g = tGrad.ref();
59 
60  const labelUList& owner = mesh.owner();
61  const labelUList& neighbour = mesh.neighbour();
62 
63  const volVectorField& C = mesh.C();
64  const surfaceVectorField& Cf = mesh.Cf();
65 
66  scalar rk = (1.0/k_ - 1.0);
67 
68  forAll(owner, facei)
69  {
70  label own = owner[facei];
71  label nei = neighbour[facei];
72 
73  scalar vsfOwn = vsf[own];
74  scalar vsfNei = vsf[nei];
75 
76  scalar maxFace = max(vsfOwn, vsfNei);
77  scalar minFace = min(vsfOwn, vsfNei);
78 
79  if (k_ < 1.0)
80  {
81  scalar maxMinFace = rk*(maxFace - minFace);
82  maxFace += maxMinFace;
83  minFace -= maxMinFace;
84  }
85 
86  // owner side
87  cellMDLimitedGrad<scalar>::limitFace
88  (
89  g[own],
90  maxFace - vsfOwn,
91  minFace - vsfOwn,
92  Cf[facei] - C[own]
93  );
94 
95  // neighbour side
96  cellMDLimitedGrad<scalar>::limitFace
97  (
98  g[nei],
99  maxFace - vsfNei,
100  minFace - vsfNei,
101  Cf[facei] - C[nei]
102  );
103  }
104 
105  const volScalarField::Boundary& bsf = vsf.boundaryField();
106 
107  forAll(bsf, patchi)
108  {
109  const fvPatchScalarField& psf = bsf[patchi];
110 
111  const labelUList& pOwner = mesh.boundary()[patchi].faceCells();
112  const vectorField& pCf = Cf.boundaryField()[patchi];
113 
114  if (psf.coupled())
115  {
116  const scalarField psfNei(psf.patchNeighbourField());
117 
118  forAll(pOwner, pFacei)
119  {
120  label own = pOwner[pFacei];
121 
122  scalar vsfOwn = vsf[own];
123  scalar vsfNei = psfNei[pFacei];
124 
125  scalar maxFace = max(vsfOwn, vsfNei);
126  scalar minFace = min(vsfOwn, vsfNei);
127 
128  if (k_ < 1.0)
129  {
130  scalar maxMinFace = rk*(maxFace - minFace);
131  maxFace += maxMinFace;
132  minFace -= maxMinFace;
133  }
134 
135  cellMDLimitedGrad<scalar>::limitFace
136  (
137  g[own],
138  maxFace - vsfOwn,
139  minFace - vsfOwn,
140  pCf[pFacei] - C[own]
141  );
142  }
143  }
144  else if (psf.fixesValue())
145  {
146  forAll(pOwner, pFacei)
147  {
148  label own = pOwner[pFacei];
149 
150  scalar vsfOwn = vsf[own];
151  scalar vsfNei = psf[pFacei];
152 
153  scalar maxFace = max(vsfOwn, vsfNei);
154  scalar minFace = min(vsfOwn, vsfNei);
155 
156  if (k_ < 1.0)
157  {
158  scalar maxMinFace = rk*(maxFace - minFace);
159  maxFace += maxMinFace;
160  minFace -= maxMinFace;
161  }
162 
163  cellMDLimitedGrad<scalar>::limitFace
164  (
165  g[own],
166  maxFace - vsfOwn,
167  minFace - vsfOwn,
168  pCf[pFacei] - C[own]
169  );
170  }
171  }
172  }
173 
174  g.correctBoundaryConditions();
176 
177  return tGrad;
178 }
179 
180 
181 template<>
184 (
185  const volVectorField& vvf,
186  const word& name
187 ) const
188 {
189  const fvMesh& mesh = vvf.mesh();
190 
191  tmp<volTensorField> tGrad = basicGradScheme_().calcGrad(vvf, name);
192 
193  if (k_ < small)
194  {
195  return tGrad;
196  }
197 
198  volTensorField& g = tGrad.ref();
199 
200  const labelUList& owner = mesh.owner();
201  const labelUList& neighbour = mesh.neighbour();
202 
203  const volVectorField& C = mesh.C();
204  const surfaceVectorField& Cf = mesh.Cf();
205 
206  scalar rk = (1.0/k_ - 1.0);
207 
208  forAll(owner, facei)
209  {
210  label own = owner[facei];
211  label nei = neighbour[facei];
212 
213  vector vvfOwn = vvf[own];
214  vector vvfNei = vvf[nei];
215 
216  vector maxFace = max(vvfOwn, vvfNei);
217  vector minFace = min(vvfOwn, vvfNei);
218 
219  if (k_ < 1.0)
220  {
221  vector maxMinFace = rk*(maxFace - minFace);
222  maxFace += maxMinFace;
223  minFace -= maxMinFace;
224  }
225 
226  // owner side
228  (
229  g[own],
230  maxFace - vvfOwn,
231  minFace - vvfOwn,
232  Cf[facei] - C[own]
233  );
234 
235 
236  // neighbour side
238  (
239  g[nei],
240  maxFace - vvfNei,
241  minFace - vvfNei,
242  Cf[facei] - C[nei]
243  );
244  }
245 
246 
247  const volVectorField::Boundary& bvf = vvf.boundaryField();
248 
249  forAll(bvf, patchi)
250  {
251  const fvPatchVectorField& psf = bvf[patchi];
252 
253  const labelUList& pOwner = mesh.boundary()[patchi].faceCells();
254  const vectorField& pCf = Cf.boundaryField()[patchi];
255 
256  if (psf.coupled())
257  {
258  const vectorField psfNei(psf.patchNeighbourField());
259 
260  forAll(pOwner, pFacei)
261  {
262  label own = pOwner[pFacei];
263 
264  vector vvfOwn = vvf[own];
265  vector vvfNei = psfNei[pFacei];
266 
267  vector maxFace = max(vvfOwn, vvfNei);
268  vector minFace = min(vvfOwn, vvfNei);
269 
270  if (k_ < 1.0)
271  {
272  vector maxMinFace = rk*(maxFace - minFace);
273  maxFace += maxMinFace;
274  minFace -= maxMinFace;
275  }
276 
278  (
279  g[own],
280  maxFace - vvfOwn, minFace - vvfOwn,
281  pCf[pFacei] - C[own]
282  );
283  }
284  }
285  else if (psf.fixesValue())
286  {
287  forAll(pOwner, pFacei)
288  {
289  label own = pOwner[pFacei];
290 
291  vector vvfOwn = vvf[own];
292  vector vvfNei = psf[pFacei];
293 
294  vector maxFace = max(vvfOwn, vvfNei);
295  vector minFace = min(vvfOwn, vvfNei);
296 
297  if (k_ < 1.0)
298  {
299  vector maxMinFace = rk*(maxFace - minFace);
300  maxFace += maxMinFace;
301  minFace -= maxMinFace;
302  }
303 
305  (
306  g[own],
307  maxFace - vvfOwn,
308  minFace - vvfOwn,
309  pCf[pFacei] - C[own]
310  );
311  }
312  }
313  }
314 
317 
318  return tGrad;
319 }
320 
321 
322 // ************************************************************************* //
static const Foam::dimensionedScalar C("C", Foam::dimTemperature, 234.5)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Graphite solid properties.
Definition: C.H:51
const Mesh & mesh() const
Return mesh.
Generic GeometricBoundaryField class.
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void correctBoundaryConditions()
Correct boundary field.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const volVectorField & C() const
Return cell centres.
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:469
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:857
const surfaceVectorField & Cf() const
Return face centres.
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:475
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
virtual bool fixesValue() const
Return true if this patch field fixes a value.
Definition: fvPatchField.H:315
virtual bool coupled() const
Return true if this patch field is coupled.
Definition: fvPatchField.H:328
virtual tmp< Field< Type > > patchNeighbourField(const Pstream::commsTypes commsType=Pstream::commsTypes::blocking) const
Return patchField on the opposite patch of a coupled patch.
Definition: fvPatchField.H:433
cellMDLimitedGrad 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.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
makeFvGradScheme(faceMDLimitedGrad)
volScalarField scalarField(fieldObject, mesh)
volVectorField vectorField(fieldObject, mesh)
label patchi
U correctBoundaryConditions()
VolField< vector > volVectorField
Definition: volFieldsFwd.H:65
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
SurfaceField< vector > surfaceVectorField
UList< label > labelUList
Definition: UList.H:65
fvPatchField< scalar > fvPatchScalarField