linearUpwindV.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 "linearUpwindV.H"
27 #include "fvMesh.H"
28 #include "volFields.H"
29 #include "surfaceFields.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 template<class Type>
36 (
38 ) const
39 {
40  const fvMesh& mesh = this->mesh();
41 
43  (
45  (
46  IOobject
47  (
48  "linearUpwindV::correction(" + vf.name() + ')',
49  mesh.time().timeName(),
50  mesh,
51  IOobject::NO_READ,
52  IOobject::NO_WRITE,
53  false
54  ),
55  mesh,
57  (
58  vf.name(),
59  vf.dimensions(),
60  Zero
61  )
62  )
63  );
64 
66 
67  const surfaceScalarField& faceFlux = this->faceFlux_;
68  const surfaceScalarField& w = mesh.weights();
69 
70  const labelList& own = mesh.owner();
71  const labelList& nei = mesh.neighbour();
72 
73  const volVectorField& C = mesh.C();
74  const surfaceVectorField& Cf = mesh.Cf();
75 
76  tmp
77  <
79  <
82  volMesh
83  >
84  > tgradVf = gradScheme_().grad(vf, gradSchemeName_);
85 
86  const GeometricField
87  <
88  typename outerProduct<vector, Type>::type,
89  fvPatchField,
90  volMesh
91  >& gradVf = tgradVf();
92 
93  forAll(faceFlux, facei)
94  {
95  vector maxCorr;
96 
97  if (faceFlux[facei] > 0.0)
98  {
99  maxCorr =
100  (1.0 - w[facei])*(vf[nei[facei]] - vf[own[facei]]);
101 
102  sfCorr[facei] =
103  (Cf[facei] - C[own[facei]]) & gradVf[own[facei]];
104  }
105  else
106  {
107  maxCorr =
108  w[facei]*(vf[own[facei]] - vf[nei[facei]]);
109 
110  sfCorr[facei] =
111  (Cf[facei] - C[nei[facei]]) & gradVf[nei[facei]];
112  }
113 
114  scalar sfCorrs = magSqr(sfCorr[facei]);
115  scalar maxCorrs = sfCorr[facei] & maxCorr;
116 
117  if (sfCorrs > 0)
118  {
119  if (maxCorrs < 0)
120  {
121  sfCorr[facei] = Zero;
122  }
123  else if (sfCorrs > maxCorrs)
124  {
125  sfCorr[facei] *= maxCorrs/(sfCorrs + vSmall);
126  }
127  }
128  }
129 
130 
132  Boundary& bSfCorr = sfCorr.boundaryFieldRef();
133 
134  forAll(bSfCorr, patchi)
135  {
136  fvsPatchField<Type>& pSfCorr = bSfCorr[patchi];
137 
138  if (pSfCorr.coupled())
139  {
140  const labelUList& pOwner =
141  mesh.boundary()[patchi].faceCells();
142 
143  const vectorField& pCf = Cf.boundaryField()[patchi];
144  const scalarField& pW = w.boundaryField()[patchi];
145 
146  const scalarField& pFaceFlux = faceFlux.boundaryField()[patchi];
147 
149  (
150  gradVf.boundaryField()[patchi].patchNeighbourField()
151  );
152 
153  const Field<Type> pVfNei
154  (
155  vf.boundaryField()[patchi].patchNeighbourField()
156  );
157 
158  // Build the d-vectors
159  vectorField pd(Cf.boundaryField()[patchi].patch().delta());
160 
161  forAll(pOwner, facei)
162  {
163  label own = pOwner[facei];
164 
165  vector maxCorr;
166 
167  if (pFaceFlux[facei] > 0)
168  {
169  pSfCorr[facei] = (pCf[facei] - C[own]) & gradVf[own];
170 
171  maxCorr = (1.0 - pW[facei])*(pVfNei[facei] - vf[own]);
172  }
173  else
174  {
175  pSfCorr[facei] =
176  (pCf[facei] - pd[facei] - C[own]) & pGradVfNei[facei];
177 
178  maxCorr = pW[facei]*(vf[own] - pVfNei[facei]);
179  }
180 
181  scalar pSfCorrs = magSqr(pSfCorr[facei]);
182  scalar maxCorrs = pSfCorr[facei] & maxCorr;
183 
184  if (pSfCorrs > 0)
185  {
186  if (maxCorrs < 0)
187  {
188  pSfCorr[facei] = Zero;
189  }
190  else if (pSfCorrs > maxCorrs)
191  {
192  pSfCorr[facei] *= maxCorrs/(pSfCorrs + vSmall);
193  }
194  }
195  }
196  }
197  }
198 
199  return tsfCorr;
200 }
201 
202 
203 namespace Foam
204 {
206 }
207 
208 
209 // ************************************************************************* //
Foam::surfaceFields.
Graphite solid properties.
Definition: C.H:48
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
const word & name() const
Return name.
Definition: IOobject.H:295
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
typeOfRank< typename pTraits< arg1 >::cmptType, direction(pTraits< arg1 >::rank)+direction(pTraits< arg2 >::rank) >::type type
Definition: products.H:90
Generic GeometricField class.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
Generic dimensioned Type class.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:284
virtual tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > correction(const GeometricField< Type, fvPatchField, volMesh > &) const
Return the explicit correction to the face-interpolate.
Definition: linearUpwindV.C:36
const dimensionSet & dimensions() const
Return dimensions.
Mesh data needed to do the Finite Volume discretisation.
Definition: volMesh.H:53
const surfaceScalarField & weights() const
Return reference to linear difference weighting factors.
dynamicFvMesh & mesh
virtual bool coupled() const
Return true if this patch field is coupled.
static const zero Zero
Definition: zero.H:97
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:278
Internal & ref()
Return a reference to the dimensioned internal field.
label patchi
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const volVectorField & C() const
Return cell centres as volVectorField.
A class for managing temporary objects.
Definition: PtrList.H:53
#define makelimitedSurfaceInterpolationTypeScheme(SS, Type)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
linearUpwindV interpolation scheme class derived from upwind and returns upwind weighting factors but...
Definition: linearUpwindV.H:51
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
Namespace for OpenFOAM.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540