FreeFOAM The Cross-Platform CFD Toolkit
boundedBackwardDdtScheme.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::fv::boundedBackwardDdtScheme
26 
27 Description
28  Second-order bounded-backward-differencing ddt using the current and
29  two previous time-step values.
30 
31 SourceFiles
32  boundedBackwardDdtScheme.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef boundedBackwardDdtScheme_H
37 #define boundedBackwardDdtScheme_H
38 
39 #include <finiteVolume/ddtScheme.H>
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace fv
50 {
51 
52 /*---------------------------------------------------------------------------*\
53  Class boundedBackwardDdtScheme Declaration
54 \*---------------------------------------------------------------------------*/
55 
57 :
58  public fv::ddtScheme<scalar>
59 {
60  // Private Member Functions
61 
62  //- Return the current time-step
63  scalar deltaT_() const;
64 
65  //- Return the previous time-step
66  scalar deltaT0_() const;
67 
68  //- Return the previous time-step or GREAT if the old timestep field
69  // wasn't available in which case Euler ddt is used
70  template<class GeoField>
71  scalar deltaT0_(const GeoField& vf) const
72  {
73  if (vf.oldTime().timeIndex() == vf.oldTime().oldTime().timeIndex())
74  {
75  return GREAT;
76  }
77  else
78  {
79  return deltaT0_();
80  }
81  }
82 
83 
84  //- Disallow default bitwise copy construct
86 
87  //- Disallow default bitwise assignment
88  void operator=(const boundedBackwardDdtScheme&);
89 
90 
91 public:
92 
93  //- Runtime type information
94  TypeName("boundedBackward");
95 
96 
97  // Constructors
98 
99  //- Construct from mesh
101  :
102  ddtScheme<scalar>(mesh)
103  {}
104 
105  //- Construct from mesh and Istream
107  :
108  ddtScheme<scalar>(mesh, is)
109  {}
110 
111 
112  // Member Functions
113 
114  //- Return mesh reference
115  const fvMesh& mesh() const
116  {
118  }
119 
121  (
122  const dimensionedScalar&
123  );
124 
126  (
127  const volScalarField&
128  );
129 
131  (
132  const dimensionedScalar&,
133  const volScalarField&
134  );
135 
137  (
138  const volScalarField&,
139  const volScalarField&
140  );
141 
143  (
145  );
146 
148  (
149  const dimensionedScalar&,
151  );
152 
154  (
155  const volScalarField&,
157  );
158 
160  (
161  const volScalarField& rA,
162  const volScalarField& U,
163  const surfaceScalarField& phi
164  );
165 
167  (
168  const volScalarField& rA,
169  const volScalarField& rho,
170  const volScalarField& U,
171  const surfaceScalarField& phi
172  );
173 
175  (
176  const volScalarField&
177  );
178 };
179 
180 
181 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
182 
183 } // End namespace fv
184 
185 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
186 
187 } // End namespace Foam
188 
189 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
190 
191 #endif
192 
193 // ************************ vim: set sw=4 sts=4 et: ************************ //