FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 
45 #include <finiteVolume/linear.H>
46 #include <finiteVolume/gaussGrad.H>
49 
50 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 
52 namespace Foam
53 {
54 
55 /*---------------------------------------------------------------------------*\
56  Class outletStabilised Declaration
57 \*---------------------------------------------------------------------------*/
58 
59 template<class Type>
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  :
95  faceFlux_
96  (
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  :
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
184  correction
185  (
187  ) const
188  {
189  if (tScheme_().corrected())
190  {
192  tScheme_().correction(vf);
193 
195 
196  const fvMesh& mesh_ = this->mesh();
197  const cellList& cells = mesh_.cells();
198 
200  {
201  if
202  (
204  (vf.boundaryField()[patchi])
206  (vf.boundaryField()[patchi])
207  )
208  {
209  const labelList& pFaceCells =
210  mesh_.boundary()[patchi].faceCells();
211 
212  forAll(pFaceCells, pFacei)
213  {
214  const cell& pFaceCell = cells[pFaceCells[pFacei]];
215 
216  forAll(pFaceCell, fi)
217  {
218  label facei = pFaceCell[fi];
219 
220  if (mesh_.isInternalFace(facei))
221  {
222  // Remove correction
223  corr[facei] = pTraits<Type>::zero;
224  }
225  }
226  }
227  }
228  }
229 
230  return tcorr;
231  }
232  else
233  {
234  return
236  }
237  }
238 };
239 
240 
241 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
242 
243 } // End namespace Foam
244 
245 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
246 
247 #endif
248 
249 // ************************ vim: set sw=4 sts=4 et: ************************ //