boundedDdtScheme.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) 2012-2016 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::fv::boundedDdtScheme
26 
27 Description
28  Bounded form of the selected ddt scheme.
29 
30  Boundedness is achieved by subtracting ddt(phi)*vf or Sp(ddt(rho), vf)
31  which is non-conservative if ddt(rho) != 0 but conservative otherwise.
32 
33  Can be used for the ddt of bounded scalar properties to improve stability
34  if insufficient convergence of the pressure equation causes temporary
35  divergence of the flux field.
36 
37 SourceFiles
38  boundedDdtScheme.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef boundedDdtScheme_H
43 #define boundedDdtScheme_H
44 
45 #include "ddtScheme.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace fv
55 {
56 
57 /*---------------------------------------------------------------------------*\
58  Class boundedDdtScheme Declaration
59 \*---------------------------------------------------------------------------*/
60 
61 template<class Type>
62 class boundedDdtScheme
63 :
64  public fv::ddtScheme<Type>
65 {
66  // Private data
67 
68  tmp<fv::ddtScheme<Type>> scheme_;
69 
70 
71  // Private Member Functions
72 
73  //- Disallow default bitwise copy construct
75 
76  //- Disallow default bitwise assignment
77  void operator=(const boundedDdtScheme&);
78 
79 
80 public:
81 
82  //- Runtime type information
83  TypeName("bounded");
84 
85 
86  // Constructors
87 
88  //- Construct from mesh and Istream
89  boundedDdtScheme(const fvMesh& mesh, Istream& is)
90  :
91  ddtScheme<Type>(mesh, is),
92  scheme_
93  (
94  fv::ddtScheme<Type>::New(mesh, is)
95  )
96  {}
97 
98 
99  // Member Functions
100 
101  //- Return mesh reference
102  const fvMesh& mesh() const
103  {
104  return fv::ddtScheme<Type>::mesh();
105  }
106 
108  (
109  const dimensioned<Type>&
110  );
111 
113  (
115  );
116 
118  (
119  const dimensionedScalar&,
121  );
122 
124  (
125  const volScalarField&,
127  );
128 
130  (
131  const volScalarField& alpha,
132  const volScalarField& rho,
134  );
135 
137  (
139  );
140 
142  (
143  const dimensionedScalar&,
145  );
146 
148  (
149  const volScalarField&,
151  );
152 
154  (
155  const volScalarField& alpha,
156  const volScalarField& rho,
158  );
161 
163  (
166  );
167 
169  (
171  const fluxFieldType& phi
172  );
173 
175  (
176  const volScalarField& rho,
179  );
180 
182  (
183  const volScalarField& rho,
185  const fluxFieldType& phi
186  );
187 
189  (
191  );
192 };
193 
194 
195 template<>
197 (
200 );
201 
202 template<>
204 (
205  const volScalarField& U,
206  const surfaceScalarField& phi
207 );
208 
209 template<>
211 (
212  const volScalarField& rho,
213  const volScalarField& U,
214  const surfaceScalarField& Uf
215 );
216 
217 template<>
219 (
220  const volScalarField& rho,
221  const volScalarField& U,
222  const surfaceScalarField& phi
223 );
224 
225 
226 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
227 
228 } // End namespace fv
229 
230 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
231 
232 } // End namespace Foam
233 
234 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
235 
236 #ifdef NoRepository
237  #include "boundedDdtScheme.C"
238 #endif
239 
240 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
241 
242 #endif
243 
244 // ************************************************************************* //
tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
surfaceScalarField & phi
Bounded form of the selected ddt scheme.
const fvMesh & mesh() const
Return mesh reference.
U
Definition: pEqn.H:83
tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
Uf
Definition: pEqn.H:67
Generic GeometricField class.
Generic dimensioned Type class.
TypeName("bounded")
Runtime type information.
Abstract base class for ddt schemes.
Definition: ddtScheme.H:64
const fvMesh & mesh() const
Return mesh reference.
Definition: ddtScheme.H:135
labelList fv(nPoints)
ddtScheme< Type >::fluxFieldType fluxFieldType
tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
static tmp< ddtScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return a pointer to a new ddtScheme created on freestore.
Definition: ddtScheme.C:46
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const volScalarField & psi
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
Namespace for OpenFOAM.
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)