FreeFOAM The Cross-Platform CFD Toolkit
ddtScheme.C
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 \*---------------------------------------------------------------------------*/
25 
26 #include <finiteVolume/fv.H>
27 #include <OpenFOAM/HashTable.H>
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace fv
38 {
39 
40 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
41 
42 template<class Type>
43 tmp<ddtScheme<Type> > ddtScheme<Type>::New
44 (
45  const fvMesh& mesh,
46  Istream& schemeData
47 )
48 {
49  if (fv::debug)
50  {
51  Info<< "ddtScheme<Type>::New(const fvMesh&, Istream&) : "
52  "constructing ddtScheme<Type>"
53  << endl;
54  }
55 
56  if (schemeData.eof())
57  {
59  (
60  "ddtScheme<Type>::New(const fvMesh&, Istream&)",
61  schemeData
62  ) << "Ddt scheme not specified" << endl << endl
63  << "Valid ddt schemes are :" << endl
64  << IstreamConstructorTablePtr_->sortedToc()
65  << exit(FatalIOError);
66  }
67 
68  word schemeName(schemeData);
69 
70  typename IstreamConstructorTable::iterator cstrIter =
71  IstreamConstructorTablePtr_->find(schemeName);
72 
73  if (cstrIter == IstreamConstructorTablePtr_->end())
74  {
76  (
77  "ddtScheme<Type>::New(const fvMesh&, Istream&)",
78  schemeData
79  ) << "unknown ddt scheme " << schemeName << endl << endl
80  << "Valid ddt schemes are :" << endl
81  << IstreamConstructorTablePtr_->sortedToc()
82  << exit(FatalIOError);
83  }
84 
85  return cstrIter()(mesh, schemeData);
86 }
87 
88 
89 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
90 
91 template<class Type>
93 {}
94 
95 
96 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
97 
98 template<class Type>
100 (
102  const fluxFieldType& phi,
103  const fluxFieldType& phiCorr
104 )
105 {
106  tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
107  - min
108  (
109  mag(phiCorr)
110  /(mag(phi) + dimensionedScalar("small", phi.dimensions(), SMALL)),
111  scalar(1)
112  );
113 
114  surfaceScalarField& ddtCouplingCoeff = tddtCouplingCoeff();
115 
117  {
118  if (U.boundaryField()[patchi].fixesValue())
119  {
120  ddtCouplingCoeff.boundaryField()[patchi] = 0.0;
121  }
122  }
123 
124  if (debug > 1)
125  {
126  Info<< "ddtCouplingCoeff mean max min = "
127  << gAverage(ddtCouplingCoeff.internalField())
128  << " " << gMax(ddtCouplingCoeff.internalField())
129  << " " << gMin(ddtCouplingCoeff.internalField())
130  << endl;
131  }
132 
133  return tddtCouplingCoeff;
134 }
135 
136 
137 template<class Type>
139 (
141  const fluxFieldType& phi
142 )
143 {
144  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
145 
146  tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
147  - min
148  (
149  mag(phi - (mesh().Sf() & fvc::interpolate(U)))
150  /(mag(phi) + dimensionedScalar("small", phi.dimensions(), VSMALL)),
151  //(rDeltaT*mesh().magSf()/mesh().deltaCoeffs()),
152  scalar(1)
153  );
154 
155  surfaceScalarField& ddtCouplingCoeff = tddtCouplingCoeff();
156 
158  {
159  if (U.boundaryField()[patchi].fixesValue())
160  {
161  ddtCouplingCoeff.boundaryField()[patchi] = 0.0;
162  }
163  }
164 
165  if (debug > 1)
166  {
167  Info<< "ddtCouplingCoeff mean max min = "
168  << gAverage(ddtCouplingCoeff.internalField())
169  << " " << gMax(ddtCouplingCoeff.internalField())
170  << " " << gMin(ddtCouplingCoeff.internalField())
171  << endl;
172  }
173 
174  return tddtCouplingCoeff;
175 }
176 
177 
178 template<class Type>
180 (
181  const volScalarField& rho,
183  const fluxFieldType& phi
184 )
185 {
186  dimensionedScalar rDeltaT = 1.0/mesh().time().deltaT();
187 
188  tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
189  - min
190  (
191  mag(phi - (mesh().Sf() & fvc::interpolate(rhoU)))
192  /(
193  mag(phi) + dimensionedScalar("small", phi.dimensions(), VSMALL)
194  //fvc::interpolate(rho)*rDeltaT
195  //*mesh().magSf()/mesh().deltaCoeffs()
196  ),
197  scalar(1)
198  );
199 
200  surfaceScalarField& ddtCouplingCoeff = tddtCouplingCoeff();
201 
202  forAll (rhoU.boundaryField(), patchi)
203  {
204  if (rhoU.boundaryField()[patchi].fixesValue())
205  {
206  ddtCouplingCoeff.boundaryField()[patchi] = 0.0;
207  }
208  }
209 
210  if (debug > 1)
211  {
212  Info<< "ddtCouplingCoeff mean max min = "
213  << gAverage(ddtCouplingCoeff.internalField())
214  << " " << gMax(ddtCouplingCoeff.internalField())
215  << " " << gMin(ddtCouplingCoeff.internalField())
216  << endl;
217  }
218 
219  return tddtCouplingCoeff;
220 }
221 
222 
223 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
224 
225 } // End namespace fv
226 
227 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
228 
229 } // End namespace Foam
230 
231 // ************************ vim: set sw=4 sts=4 et: ************************ //