ddtScheme.C
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-2024 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 \*---------------------------------------------------------------------------*/
25 
26 #include "fv.H"
27 #include "HashTable.H"
28 #include "surfaceInterpolate.H"
29 #include "fvMatrix.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 
36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37 
38 namespace fv
39 {
40 
41 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
42 
43 template<class Type>
45 (
46  const fvMesh& mesh,
47  Istream& schemeData
48 )
49 {
50  if (fv::debug)
51  {
52  InfoInFunction << "Constructing ddtScheme<Type>" << endl;
53  }
54 
55  if (schemeData.eof())
56  {
58  (
59  schemeData
60  ) << "Ddt scheme not specified" << endl << endl
61  << "Valid ddt schemes are :" << endl
62  << IstreamConstructorTablePtr_->sortedToc()
63  << exit(FatalIOError);
64  }
65 
66  const word schemeName(schemeData);
67 
68  typename IstreamConstructorTable::iterator cstrIter =
69  IstreamConstructorTablePtr_->find(schemeName);
70 
71  if (cstrIter == IstreamConstructorTablePtr_->end())
72  {
74  (
75  schemeData
76  ) << "Unknown ddt scheme " << schemeName << nl << nl
77  << "Valid ddt schemes are :" << endl
78  << IstreamConstructorTablePtr_->sortedToc()
79  << exit(FatalIOError);
80  }
81 
82  return cstrIter()(mesh, schemeData);
83 }
84 
85 
86 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
87 
88 template<class Type>
90 {}
91 
92 
93 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
94 
95 template<class Type>
97 (
98  const volScalarField& alpha,
99  const volScalarField& rho,
100  const VolField<Type>& vf
101 )
102 {
104 
105  return tmp<VolField<Type>>
106  (
108  );
109 }
110 
111 
112 template<class Type>
114 (
115  const volScalarField& alpha,
116  const volScalarField& rho,
117  const VolField<Type>& vf
118 )
119 {
121 
122  return tmp<fvMatrix<Type>>
123  (
124  new fvMatrix<Type>
125  (
126  vf,
127  alpha.dimensions()*rho.dimensions()
129  )
130  );
131 }
132 
133 
134 template<class Type>
136 (
137  const SurfaceField<Type>& sf
138 )
139 {
141 
142  return tmp<SurfaceField<Type>>
143  (
145  );
146 }
147 
148 
149 
150 template<class Type>
152 (
153  const VolField<Type>& U,
154  const fluxFieldType& phi,
155  const fluxFieldType& phiCorr
156 )
157 {
158  // Courant number limited formulation
159  /*
160  tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
161  - min
162  (
163  mag(phiCorr)
164  *mesh().time().deltaT()*mesh().deltaCoeffs()/mesh().magSf(),
165  scalar(1)
166  );
167  */
168 
169  // Flux normalised formulation
170  tmp<surfaceScalarField> tddtCouplingCoeff = scalar(1)
171  - min
172  (
173  mag(phiCorr)
174  /(mag(phi) + dimensionedScalar("small", phi.dimensions(), SMALL)),
175  scalar(1)
176  );
177 
178  surfaceScalarField& ddtCouplingCoeff = tddtCouplingCoeff.ref();
179 
180  surfaceScalarField::Boundary& ccbf = ddtCouplingCoeff.boundaryFieldRef();
181 
182  forAll(U.boundaryField(), patchi)
183  {
184  if (!U.boundaryField()[patchi].coupled())
185  {
186  ccbf[patchi] = 0;
187  }
188  }
189 
190  if (debug > 1)
191  {
193  << "ddtCouplingCoeff mean max min = "
194  << gAverage(ddtCouplingCoeff.primitiveField())
195  << " " << gMax(ddtCouplingCoeff.primitiveField())
196  << " " << gMin(ddtCouplingCoeff.primitiveField())
197  << endl;
198  }
199 
200  return tddtCouplingCoeff;
201 }
202 
203 
204 template<class Type>
206 (
207  const VolField<Type>& U,
208  const fluxFieldType& phi,
209  const fluxFieldType& phiCorr,
210  const volScalarField& rho
211 )
212 {
213  return fvcDdtPhiCoeff(U, phi, phiCorr);
214 }
215 
216 
217 template<class Type>
219 (
220  const VolField<Type>& U,
221  const fluxFieldType& phi
222 )
223 {
224  return fvcDdtPhiCoeff(U, phi, phi - fvc::dotInterpolate(mesh().Sf(), U));
225 }
226 
227 
228 template<class Type>
230 (
231  const VolField<Type>& U,
232  const fluxFieldType& phi,
233  const volScalarField& rho
234 )
235 {
236  return fvcDdtPhiCoeff
237  (
238  U,
239  phi,
240  (phi - fvc::dotInterpolate(mesh().Sf(), U))
241  );
242 }
243 
244 
245 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
246 
247 } // End namespace fv
248 
249 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 
251 } // End namespace Foam
252 
253 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const dimensionSet & dimensions() const
Return dimensions.
Generic GeometricBoundaryField class.
Generic GeometricField class.
const Internal::FieldType & primitiveField() const
Return a const-reference to the primitive field.
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
bool eof() const
Return true if end of input seen.
Definition: IOstream.H:336
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
virtual tmp< fvMatrix< Type > > fvmDdt(const VolField< Type > &)=0
virtual ~ddtScheme()
Destructor.
Definition: ddtScheme.C:89
virtual tmp< VolField< Type > > fvcDdt(const dimensioned< Type > &)=0
virtual tmp< surfaceScalarField > fvcDdtPhiCoeff(const VolField< Type > &U, const fluxFieldType &phi, const fluxFieldType &phiCorr)
Definition: ddtScheme.C:152
static tmp< ddtScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return a pointer to a new ddtScheme created on freestore.
Definition: ddtScheme.C:45
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:346
label patchi
volScalarField sf(fieldObject, mesh)
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
#define InfoInFunction
Report an information message using Foam::Info.
static tmp< SurfaceField< typename innerProduct< vector, Type >::type > > dotInterpolate(const surfaceVectorField &Sf, const VolField< Type > &tvf)
Interpolate field onto faces.
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static const scalar SMALL
Definition: scalar.H:115
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimTime
const dimensionSet dimVolume
dimensioned< scalar > mag(const dimensioned< Type > &)
IOerror FatalIOError
Type gAverage(const FieldField< Field, Type > &f)
Type gMin(const FieldField< Field, Type > &f)
static const char nl
Definition: Ostream.H:266
Type gMax(const FieldField< Field, Type > &f)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
labelList fv(nPoints)