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-2019 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  "linearUpwindV::correction(" + vf.name() + ')',
47  mesh,
49  (
50  vf.name(),
51  vf.dimensions(),
52  Zero
53  )
54  )
55  );
56 
58 
59  const surfaceScalarField& faceFlux = this->faceFlux_;
60  const surfaceScalarField& w = mesh.weights();
61 
62  const labelList& own = mesh.owner();
63  const labelList& nei = mesh.neighbour();
64 
65  const volVectorField& C = mesh.C();
66  const surfaceVectorField& Cf = mesh.Cf();
67 
68  tmp
69  <
71  <
74  volMesh
75  >
76  > tgradVf = gradScheme_().grad(vf, gradSchemeName_);
77 
78  const GeometricField
79  <
80  typename outerProduct<vector, Type>::type,
81  fvPatchField,
82  volMesh
83  >& gradVf = tgradVf();
84 
85  forAll(faceFlux, facei)
86  {
87  vector maxCorr;
88 
89  if (faceFlux[facei] > 0.0)
90  {
91  maxCorr =
92  (1.0 - w[facei])*(vf[nei[facei]] - vf[own[facei]]);
93 
94  sfCorr[facei] =
95  (Cf[facei] - C[own[facei]]) & gradVf[own[facei]];
96  }
97  else
98  {
99  maxCorr =
100  w[facei]*(vf[own[facei]] - vf[nei[facei]]);
101 
102  sfCorr[facei] =
103  (Cf[facei] - C[nei[facei]]) & gradVf[nei[facei]];
104  }
105 
106  scalar sfCorrs = magSqr(sfCorr[facei]);
107  scalar maxCorrs = sfCorr[facei] & maxCorr;
108 
109  if (sfCorrs > 0)
110  {
111  if (maxCorrs < 0)
112  {
113  sfCorr[facei] = Zero;
114  }
115  else if (sfCorrs > maxCorrs)
116  {
117  sfCorr[facei] *= maxCorrs/(sfCorrs + vSmall);
118  }
119  }
120  }
121 
122 
124  Boundary& bSfCorr = sfCorr.boundaryFieldRef();
125 
126  forAll(bSfCorr, patchi)
127  {
128  fvsPatchField<Type>& pSfCorr = bSfCorr[patchi];
129 
130  if (pSfCorr.coupled())
131  {
132  const labelUList& pOwner =
133  mesh.boundary()[patchi].faceCells();
134 
135  const vectorField& pCf = Cf.boundaryField()[patchi];
136  const scalarField& pW = w.boundaryField()[patchi];
137 
138  const scalarField& pFaceFlux = faceFlux.boundaryField()[patchi];
139 
141  (
142  gradVf.boundaryField()[patchi].patchNeighbourField()
143  );
144 
145  const Field<Type> pVfNei
146  (
147  vf.boundaryField()[patchi].patchNeighbourField()
148  );
149 
150  // Build the d-vectors
151  vectorField pd(Cf.boundaryField()[patchi].patch().delta());
152 
153  forAll(pOwner, facei)
154  {
155  label own = pOwner[facei];
156 
157  vector maxCorr;
158 
159  if (pFaceFlux[facei] > 0)
160  {
161  pSfCorr[facei] = (pCf[facei] - C[own]) & gradVf[own];
162 
163  maxCorr = (1.0 - pW[facei])*(pVfNei[facei] - vf[own]);
164  }
165  else
166  {
167  pSfCorr[facei] =
168  (pCf[facei] - pd[facei] - C[own]) & pGradVfNei[facei];
169 
170  maxCorr = pW[facei]*(vf[own] - pVfNei[facei]);
171  }
172 
173  scalar pSfCorrs = magSqr(pSfCorr[facei]);
174  scalar maxCorrs = pSfCorr[facei] & maxCorr;
175 
176  if (pSfCorrs > 0)
177  {
178  if (maxCorrs < 0)
179  {
180  pSfCorr[facei] = Zero;
181  }
182  else if (pSfCorrs > maxCorrs)
183  {
184  pSfCorr[facei] *= maxCorrs/(pSfCorrs + vSmall);
185  }
186  }
187  }
188  }
189  }
190 
191  return tsfCorr;
192 }
193 
194 
195 namespace Foam
196 {
198 }
199 
200 
201 // ************************************************************************* //
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:303
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.
Generic dimensioned Type class.
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)
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