outletStabilised.H
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-2020 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"
48 #include "mixedFvPatchField.H"
50 
51 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
52 
53 namespace Foam
54 {
55 
56 /*---------------------------------------------------------------------------*\
57  Class outletStabilised Declaration
58 \*---------------------------------------------------------------------------*/
59 
60 template<class Type>
61 class outletStabilised
62 :
63  public surfaceInterpolationScheme<Type>
64 {
65  // Private member data
66 
67  const surfaceScalarField& faceFlux_;
69 
70 
71 public:
72 
73  //- Runtime type information
74  TypeName("outletStabilised");
75 
76 
77  // Constructors
78 
79  //- Construct from mesh and Istream
81  (
82  const fvMesh& mesh,
83  Istream& is
84  )
85  :
86  surfaceInterpolationScheme<Type>(mesh),
87  faceFlux_
88  (
89  mesh.lookupObject<surfaceScalarField>
90  (
91  word(is)
92  )
93  ),
94  tScheme_
95  (
96  surfaceInterpolationScheme<Type>::New(mesh, faceFlux_, is)
97  )
98  {}
99 
100 
101  //- Construct from mesh, faceFlux and Istream
103  (
104  const fvMesh& mesh,
105  const surfaceScalarField& faceFlux,
106  Istream& is
107  )
108  :
109  surfaceInterpolationScheme<Type>(mesh),
110  faceFlux_(faceFlux),
111  tScheme_
112  (
113  surfaceInterpolationScheme<Type>::New(mesh, faceFlux, is)
114  )
115  {}
116 
117  //- Disallow default bitwise copy construction
118  outletStabilised(const outletStabilised&) = delete;
119 
120 
121  // Member Functions
122 
123  //- Return the interpolation weighting factors
125  (
127  ) const
128  {
129  tmp<surfaceScalarField> tw = tScheme_().weights(vf);
130  surfaceScalarField& w = tw.ref();
131 
132  const fvMesh& mesh_ = this->mesh();
133  const cellList& cells = mesh_.cells();
134 
136  {
137  if
138  (
140  (vf.boundaryField()[patchi])
143  (vf.boundaryField()[patchi])
144  )
145  {
146  const labelList& pFaceCells =
147  mesh_.boundary()[patchi].faceCells();
148 
149  forAll(pFaceCells, pFacei)
150  {
151  const cell& pFaceCell = cells[pFaceCells[pFacei]];
152 
153  forAll(pFaceCell, fi)
154  {
155  label facei = pFaceCell[fi];
156 
157  if (mesh_.isInternalFace(facei))
158  {
159  // Apply upwind differencing
160  w[facei] = pos0(faceFlux_[facei]);
161  }
162  }
163  }
164  }
165  }
166 
167  return tw;
168  }
169 
170  //- Return true if this scheme uses an explicit correction
171  virtual bool corrected() const
172  {
173  return tScheme_().corrected();
174  }
175 
176  //- Return the explicit correction to the face-interpolate
177  // set to zero on the near-boundary faces where upwind is applied
180  (
182  ) const
183  {
184  if (tScheme_().corrected())
185  {
187  tScheme_().correction(vf);
188 
190  tcorr.ref();
191 
192  const fvMesh& mesh_ = this->mesh();
193  const cellList& cells = mesh_.cells();
194 
196  {
197  if
198  (
200  (vf.boundaryField()[patchi])
202  (vf.boundaryField()[patchi])
203  )
204  {
205  const labelList& pFaceCells =
206  mesh_.boundary()[patchi].faceCells();
207 
208  forAll(pFaceCells, pFacei)
209  {
210  const cell& pFaceCell = cells[pFaceCells[pFacei]];
211 
212  forAll(pFaceCell, fi)
213  {
214  label facei = pFaceCell[fi];
215 
216  if (mesh_.isInternalFace(facei))
217  {
218  // Remove correction
219  corr[facei] = Zero;
220  }
221  }
222  }
223  }
224  }
225 
226  return tcorr;
227  }
228  else
229  {
231  (
232  nullptr
233  );
234  }
235  }
236 
237 
238  // Member Operators
239 
240  //- Disallow default bitwise assignment
241  void operator=(const outletStabilised&) = delete;
242 };
243 
244 
245 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
246 
247 } // End namespace Foam
248 
249 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 
251 #endif
252 
253 // ************************************************************************* //
virtual bool corrected() const
Return true if this scheme uses an explicit correction.
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
#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
void operator=(const outletStabilised &)=delete
Disallow default bitwise assignment.
outletStabilised(const fvMesh &mesh, Istream &is)
Construct from mesh and Istream.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
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.
const cellList & cells() const
const fvMesh & mesh() const
Return mesh reference.
TypeName("outletStabilised")
Runtime type information.
tmp< surfaceScalarField > weights(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the interpolation weighting factors.
Base class for direction-mixed boundary conditions.
Outlet-stabilised interpolation scheme which applies upwind differencing to the faces of the cells ad...
const cellShapeList & cells
A class for handling words, derived from string.
Definition: word.H:59
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
static const zero Zero
Definition: zero.H:97
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
dimensionedScalar pos0(const dimensionedScalar &ds)
This boundary condition applies a zero-gradient condition from the patch internal field onto the patc...
label patchi
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A class for managing temporary objects.
Definition: PtrList.H:53
Abstract base class for surface interpolation schemes.
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.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540