linearUpwindVTemplates.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-2022 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 "linearUpwindV.H"
27 
28 // * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
29 
30 template<class Type>
32 (
33  const fvMesh& mesh,
34  const surfaceScalarField& faceFlux
35 )
36 :
37  upwind<Type>(mesh, faceFlux),
38  gradSchemeName_("grad"),
39  gradScheme_
40  (
41  new fv::gaussGrad<Type>(mesh)
42  )
43 {}
44 
45 
46 template<class Type>
48 (
49  const fvMesh& mesh,
50  Istream& schemeData
51 )
52 :
53  upwind<Type>(mesh, schemeData),
54  gradSchemeName_(schemeData),
55  gradScheme_
56  (
57  fv::gradScheme<Type>::New
58  (
59  mesh,
60  mesh.schemes().grad(gradSchemeName_)
61  )
62  )
63 {}
64 
65 
66 template<class Type>
68 (
69  const fvMesh& mesh,
70  const surfaceScalarField& faceFlux,
71  Istream& schemeData
72 )
73 :
74  upwind<Type>(mesh, faceFlux, schemeData),
75  gradSchemeName_(schemeData),
76  gradScheme_
77  (
78  fv::gradScheme<Type>::New
79  (
80  mesh,
81  mesh.schemes().grad(gradSchemeName_)
82  )
83  )
84 {}
85 
86 
87 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
88 
89 template<class Type>
92 (
93  const VolField<Type>& vf
94 ) const
95 {
96  const fvMesh& mesh = this->mesh();
97 
98  tmp<SurfaceField<Type>> tsfCorr
99  (
101  (
102  "linearUpwindV::correction(" + vf.name() + ')',
103  mesh,
105  (
106  vf.name(),
107  vf.dimensions(),
108  Zero
109  )
110  )
111  );
112 
113  SurfaceField<Type>& sfCorr = tsfCorr.ref();
114 
115  const surfaceScalarField& faceFlux = this->faceFlux_;
116  const surfaceScalarField& w = mesh.weights();
117 
118  const labelList& own = mesh.owner();
119  const labelList& nei = mesh.neighbour();
120 
121  const volVectorField& C = mesh.C();
122  const surfaceVectorField& Cf = mesh.Cf();
123 
125  (
126  gradScheme_().grad(vf, gradSchemeName_)
127  );
128 
130  = tgradVf();
131 
132  forAll(faceFlux, facei)
133  {
134  vector maxCorr;
135 
136  if (faceFlux[facei] > 0.0)
137  {
138  maxCorr =
139  (1.0 - w[facei])*(vf[nei[facei]] - vf[own[facei]]);
140 
141  sfCorr[facei] =
142  (Cf[facei] - C[own[facei]]) & gradVf[own[facei]];
143  }
144  else
145  {
146  maxCorr =
147  w[facei]*(vf[own[facei]] - vf[nei[facei]]);
148 
149  sfCorr[facei] =
150  (Cf[facei] - C[nei[facei]]) & gradVf[nei[facei]];
151  }
152 
153  scalar sfCorrs = magSqr(sfCorr[facei]);
154  scalar maxCorrs = sfCorr[facei] & maxCorr;
155 
156  if (sfCorrs > 0)
157  {
158  if (maxCorrs < 0)
159  {
160  sfCorr[facei] = Zero;
161  }
162  else if (sfCorrs > maxCorrs)
163  {
164  sfCorr[facei] *= maxCorrs/(sfCorrs + vSmall);
165  }
166  }
167  }
168 
169 
170  typename SurfaceField<Type>::
171  Boundary& bSfCorr = sfCorr.boundaryFieldRef();
172 
173  forAll(bSfCorr, patchi)
174  {
175  fvsPatchField<Type>& pSfCorr = bSfCorr[patchi];
176 
177  if (pSfCorr.coupled())
178  {
179  const labelUList& pOwner =
180  mesh.boundary()[patchi].faceCells();
181 
182  const vectorField& pCf = Cf.boundaryField()[patchi];
183  const scalarField& pW = w.boundaryField()[patchi];
184 
185  const scalarField& pFaceFlux = faceFlux.boundaryField()[patchi];
186 
188  (
189  gradVf.boundaryField()[patchi].patchNeighbourField()
190  );
191 
192  const Field<Type> pVfNei
193  (
194  vf.boundaryField()[patchi].patchNeighbourField()
195  );
196 
197  // Build the d-vectors
198  vectorField pd(Cf.boundaryField()[patchi].patch().delta());
199 
200  forAll(pOwner, facei)
201  {
202  label own = pOwner[facei];
203 
204  vector maxCorr;
205 
206  if (pFaceFlux[facei] > 0)
207  {
208  pSfCorr[facei] = (pCf[facei] - C[own]) & gradVf[own];
209 
210  maxCorr = (1.0 - pW[facei])*(pVfNei[facei] - vf[own]);
211  }
212  else
213  {
214  pSfCorr[facei] =
215  (pCf[facei] - pd[facei] - C[own]) & pGradVfNei[facei];
216 
217  maxCorr = pW[facei]*(vf[own] - pVfNei[facei]);
218  }
219 
220  scalar pSfCorrs = magSqr(pSfCorr[facei]);
221  scalar maxCorrs = pSfCorr[facei] & maxCorr;
222 
223  if (pSfCorrs > 0)
224  {
225  if (maxCorrs < 0)
226  {
227  pSfCorr[facei] = Zero;
228  }
229  else if (pSfCorrs > maxCorrs)
230  {
231  pSfCorr[facei] *= maxCorrs/(pSfCorrs + vSmall);
232  }
233  }
234  }
235  }
236  }
237 
238  return tsfCorr;
239 }
240 
241 
242 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Graphite solid properties.
Definition: C.H:51
const dimensionSet & dimensions() const
Return dimensions.
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const word & name() const
Return name.
Definition: IOobject.H:310
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
Generic dimensioned Type class.
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
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:83
virtual bool coupled() const
Return true if this patch field is coupled.
linearUpwindV(const fvMesh &mesh, const surfaceScalarField &faceFlux)
Construct from faceFlux.
virtual tmp< SurfaceField< Type > > correction(const VolField< Type > &) const
Return the explicit correction to the face-interpolate.
const surfaceScalarField & weights() const
Return reference to linear difference weighting factors.
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
Upwind interpolation scheme class.
Definition: upwind.H:54
label patchi
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
static const zero Zero
Definition: zero.H:97
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
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
labelList fv(nPoints)