All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
age.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) 2018-2020 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 "age.H"
27 #include "fvmDiv.H"
28 #include "fvmLaplacian.H"
29 #include "momentumTransportModel.H"
31 #include "wallFvPatch.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace functionObjects
40 {
41  defineTypeNameAndDebug(age, 0);
42  addToRunTimeSelectionTable(functionObject, age, dictionary);
43 }
44 }
45 
46 
47 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
48 
49 Foam::wordList Foam::functionObjects::age::patchTypes() const
50 {
51  wordList result
52  (
53  mesh_.boundary().size(),
55  );
56 
57  forAll(mesh_.boundary(), patchi)
58  {
59  if (isA<wallFvPatch>(mesh_.boundary()[patchi]))
60  {
62  }
63  }
64 
65  return result;
66 }
67 
68 
69 bool Foam::functionObjects::age::converged
70 (
71  const int nCorr,
72  const scalar initialResidual
73 ) const
74 {
75  if (initialResidual < tolerance_)
76  {
77  Info<< "Field " << typeName
78  << " converged in " << nCorr << " correctors\n" << endl;
79 
80  return true;
81  }
82  else
83  {
84  return false;
85  }
86 }
87 
88 
89 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
90 
92 (
93  const word& name,
94  const Time& runTime,
95  const dictionary& dict
96 )
97 :
98  fvMeshFunctionObject(name, runTime, dict),
99  fvOptions_(mesh_)
100 {
101  read(dict);
102 }
103 
104 
105 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
106 
108 {}
109 
110 
111 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
112 
114 {
115  phiName_ = dict.lookupOrDefault<word>("phi", "phi");
116  rhoName_ = dict.lookupOrDefault<word>("rho", "rho");
117  nCorr_ = dict.lookupOrDefault<int>("nCorr", 5);
118  schemesField_ = dict.lookupOrDefault<word>("schemesField", typeName);
119  diffusion_ = dict.lookupOrDefault<Switch>("diffusion", false);
120  tolerance_ = dict.lookupOrDefault<scalar>("tolerance", 1e-5);
121 
122  if (dict.found("fvOptions"))
123  {
124  fvOptions_.reset(dict.subDict("fvOptions"));
125  }
126 
127  return true;
128 }
129 
130 
132 {
134  (
135  new volScalarField
136  (
137  IOobject
138  (
139  typeName,
140  mesh_.time().timeName(),
141  mesh_,
144  false
145  ),
146  mesh_,
148  patchTypes()
149  )
150  );
151  volScalarField& age = tage.ref();
152 
153  const word divScheme("div(phi," + schemesField_ + ")");
154 
155  // Set under-relaxation coeff
156  scalar relaxCoeff = 0.0;
157  if (mesh_.relaxEquation(schemesField_))
158  {
159  relaxCoeff = mesh_.equationRelaxationFactor(schemesField_);
160  }
161 
162  // This only works because the null constructed inletValue for an
163  // inletOutletFvPatchField is zero. If we needed any other value we would
164  // have to loop over the inletOutlet patches and explicitly set the
165  // inletValues. We would need to change the interface of inletOutlet in
166  // order to do this.
167 
168  const surfaceScalarField& phi =
169  mesh_.lookupObject<surfaceScalarField>(phiName_);
170 
171  if (phi.dimensions() == dimMass/dimTime)
172  {
173  const volScalarField& rho =
174  mesh_.lookupObject<volScalarField>(rhoName_);
175 
176  tmp<volScalarField> tmuEff;
177  word laplacianScheme;
178 
179  if (diffusion_)
180  {
181  tmuEff =
182  mesh_.lookupObject<momentumTransportModel>
183  (
184  momentumTransportModel::typeName
185  ).muEff();
186 
187  laplacianScheme =
188  "laplacian(" + tmuEff().name() + ',' + schemesField_ + ")";
189  }
190 
191  for (int i=0; i<=nCorr_; i++)
192  {
193  fvScalarMatrix ageEqn
194  (
195  fvm::div(phi, age, divScheme) == rho + fvOptions_(rho, age)
196  );
197 
198  if (diffusion_)
199  {
200  ageEqn -= fvm::laplacian(tmuEff(), age, laplacianScheme);
201  }
202 
203  ageEqn.relax(relaxCoeff);
204 
205  fvOptions_.constrain(ageEqn);
206 
207  if (converged(i, ageEqn.solve(schemesField_).initialResidual()))
208  {
209  break;
210  };
211  }
212  }
213  else
214  {
215  tmp<volScalarField> tnuEff;
216  word laplacianScheme;
217 
218  if (diffusion_)
219  {
220  tnuEff =
221  mesh_.lookupObject<momentumTransportModel>
222  (
223  momentumTransportModel::typeName
224  ).nuEff();
225 
226  laplacianScheme =
227  "laplacian(" + tnuEff().name() + ',' + schemesField_ + ")";
228  }
229 
230  for (int i=0; i<=nCorr_; i++)
231  {
232  fvScalarMatrix ageEqn
233  (
234  fvm::div(phi, age, divScheme)
235  == dimensionedScalar(1) + fvOptions_(age)
236  );
237 
238  if (diffusion_)
239  {
240  ageEqn -= fvm::laplacian(tnuEff(), age, laplacianScheme);
241  }
242 
243  ageEqn.relax(relaxCoeff);
244 
245  fvOptions_.constrain(ageEqn);
246 
247  if (converged(i, ageEqn.solve(schemesField_).initialResidual()))
248  {
249  break;
250  }
251  }
252  }
253 
254  Info<< "Min/max age:" << min(age).value()
255  << ' ' << max(age).value() << endl;
256 
257  return store(tage);
258 }
259 
260 
262 {
263  return writeObject(typeName);
264 }
265 
266 
267 // ************************************************************************* //
virtual bool execute()
Execute.
Definition: age.C:131
virtual bool read(const dictionary &)
Read the data.
Definition: age.C:113
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
static const char *const typeName
Definition: Field.H:105
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
age(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: age.C:92
addToRunTimeSelectionTable(functionObject, Qdot, dictionary)
wordList patchTypes(nPatches)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
Calculate the matrix for the laplacian of the field.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
Macros for easy insertion into run-time selection tables.
phi
Definition: pEqn.H:104
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:934
const dimensionSet & dimensions() const
Return dimensions.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
A class for handling words, derived from string.
Definition: word.H:59
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:521
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
defineTypeNameAndDebug(Qdot, 0)
Calculates and writes out the time taken for a particle to travel from an inlet to the location...
Definition: age.H:130
Abstract base class for turbulence models (RAS, LES and laminar).
Calculate the matrix for the divergence of the given field and flux.
List< word > wordList
A List of words.
Definition: fileName.H:54
label patchi
virtual ~age()
Destructor.
Definition: age.C:107
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
virtual bool write()
Write.
Definition: age.C:261
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.