outletStabilised.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011 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 Class
25  Foam::outletStabilised
26 
27 Description
28  Outlet-stabilised interpolation scheme which applies upwind differencing
29  to the faces of the cells adjacent to outlets.
30 
31  This is particularly useful to stabilise the velocity at entrainment
32  boundaries for LES cases using linear or other centred differencing
33  schemes.
34 
35 SourceFiles
36  outletStabilised.C
37 
38 \*---------------------------------------------------------------------------*/
39 
40 #ifndef outletStabilised_H
41 #define outletStabilised_H
42 
44 #include "skewCorrectionVectors.H"
45 #include "linear.H"
46 #include "gaussGrad.H"
47 #include "mixedFvPatchField.H"
49 
50 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 
52 namespace Foam
53 {
54 
55 /*---------------------------------------------------------------------------*\
56  Class outletStabilised Declaration
57 \*---------------------------------------------------------------------------*/
58 
59 template<class Type>
60 class outletStabilised
61 :
62  public surfaceInterpolationScheme<Type>
63 {
64  // Private member data
65 
66  const surfaceScalarField& faceFlux_;
68 
69 
70  // Private Member Functions
71 
72  //- Disallow default bitwise copy construct
74 
75  //- Disallow default bitwise assignment
76  void operator=(const outletStabilised&);
77 
78 
79 public:
80 
81  //- Runtime type information
82  TypeName("outletStabilised");
83 
84 
85  // Constructors
86 
87  //- Construct from mesh and Istream
89  (
90  const fvMesh& mesh,
91  Istream& is
92  )
93  :
94  surfaceInterpolationScheme<Type>(mesh),
95  faceFlux_
96  (
97  mesh.lookupObject<surfaceScalarField>
98  (
99  word(is)
100  )
101  ),
102  tScheme_
103  (
104  surfaceInterpolationScheme<Type>::New(mesh, faceFlux_, is)
105  )
106  {}
107 
108 
109  //- Construct from mesh, faceFlux and Istream
111  (
112  const fvMesh& mesh,
113  const surfaceScalarField& faceFlux,
114  Istream& is
115  )
116  :
117  surfaceInterpolationScheme<Type>(mesh),
118  faceFlux_(faceFlux),
119  tScheme_
120  (
121  surfaceInterpolationScheme<Type>::New(mesh, faceFlux, is)
122  )
123  {}
124 
125 
126  // Member Functions
127 
128  //- Return the interpolation weighting factors
130  (
132  ) const
133  {
134  tmp<surfaceScalarField> tw = tScheme_().weights(vf);
135  surfaceScalarField& w = tw();
136 
137  const fvMesh& mesh_ = this->mesh();
138  const cellList& cells = mesh_.cells();
139 
141  {
142  if
143  (
145  (vf.boundaryField()[patchi])
148  (vf.boundaryField()[patchi])
149  )
150  {
151  const labelList& pFaceCells =
152  mesh_.boundary()[patchi].faceCells();
153 
154  forAll(pFaceCells, pFacei)
155  {
156  const cell& pFaceCell = cells[pFaceCells[pFacei]];
157 
158  forAll(pFaceCell, fi)
159  {
160  label facei = pFaceCell[fi];
161 
162  if (mesh_.isInternalFace(facei))
163  {
164  // Apply upwind differencing
165  w[facei] = pos(faceFlux_[facei]);
166  }
167  }
168  }
169  }
170  }
171 
172  return tw;
173  }
174 
175  //- Return true if this scheme uses an explicit correction
176  virtual bool corrected() const
177  {
178  return tScheme_().corrected();
179  }
180 
181  //- Return the explicit correction to the face-interpolate
182  // set to zero on the near-boundary faces where upwinf is applied
185  (
187  ) const
188  {
189  if (tScheme_().corrected())
190  {
192  tScheme_().correction(vf);
193 
195  tcorr();
196 
197  const fvMesh& mesh_ = this->mesh();
198  const cellList& cells = mesh_.cells();
199 
201  {
202  if
203  (
205  (vf.boundaryField()[patchi])
207  (vf.boundaryField()[patchi])
208  )
209  {
210  const labelList& pFaceCells =
211  mesh_.boundary()[patchi].faceCells();
212 
213  forAll(pFaceCells, pFacei)
214  {
215  const cell& pFaceCell = cells[pFaceCells[pFacei]];
216 
217  forAll(pFaceCell, fi)
218  {
219  label facei = pFaceCell[fi];
220 
221  if (mesh_.isInternalFace(facei))
222  {
223  // Remove correction
224  corr[facei] = pTraits<Type>::zero;
225  }
226  }
227  }
228  }
229  }
230 
231  return tcorr;
232  }
233  else
234  {
236  (
237  NULL
238  );
239  }
240  }
241 };
242 
243 
244 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
245 
246 } // End namespace Foam
247 
248 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
249 
250 #endif
251 
252 // ************************************************************************* //
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
This boundary condition provides a base class for &#39;mixed&#39; type boundary conditions, i.e. conditions that mix fixed value and patch-normal gradient conditions.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
A class for handling words, derived from string.
Definition: word.H:59
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 cellList & cells() const
TypeName("outletStabilised")
Runtime type information.
virtual tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > correction(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the explicit correction to the face-interpolate.
Namespace for OpenFOAM.
This boundary condition applies a zero-gradient condition from the patch internal field onto the patc...
Abstract base class for surface interpolation schemes.
const fvMesh & mesh() const
Return mesh reference.
#define forAll(list, i)
Definition: UList.H:421
label patchi
Base class for direction-mixed boundary conditions.
const cellShapeList & cells
Traits class for primitives.
Definition: pTraits.H:50
Outlet-stabilised interpolation scheme which applies upwind differencing to the faces of the cells ad...
dimensionedScalar pos(const dimensionedScalar &ds)
tmp< surfaceScalarField > weights(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the interpolation weighting factors.
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
A class for managing temporary objects.
Definition: PtrList.H:118
virtual bool corrected() const
Return true if this scheme uses an explicit correction.
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:552