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-2023 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 "fvModels.H"
30 #include "fvConstraints.H"
31 #include "momentumTransportModel.H"
33 #include "wallFvPatch.H"
36 
37 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41 namespace functionObjects
42 {
45 }
46 }
47 
48 
49 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
51 Foam::wordList Foam::functionObjects::age::patchTypes() const
52 {
53  wordList result
54  (
55  mesh_.boundary().size(),
57  );
58 
60  {
61  if (isA<wallFvPatch>(mesh_.boundary()[patchi]))
62  {
64  }
65  }
66 
67  return result;
68 }
69 
70 
71 bool Foam::functionObjects::age::converged
72 (
73  const int nCorr,
74  const scalar initialResidual
75 ) const
76 {
77  if (initialResidual < tolerance_)
78  {
79  Info<< "Field " << typeName
80  << " converged in " << nCorr << " correctors\n" << endl;
81 
82  return true;
83  }
84  else
85  {
86  return false;
87  }
88 }
89 
90 
91 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
92 
94 (
95  const word& name,
96  const Time& runTime,
97  const dictionary& dict
98 )
99 :
100  fvMeshFunctionObject(name, runTime, dict)
101 {
102  read(dict);
103 }
104 
105 
106 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
107 
109 {}
110 
111 
112 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
113 
115 {
116  phiName_ = dict.lookupOrDefault<word>("phi", "phi");
117  rhoName_ = dict.lookupOrDefault<word>("rho", "rho");
118  nCorr_ = dict.lookupOrDefault<int>("nCorr", 5);
119  schemesField_ = dict.lookupOrDefault<word>("schemesField", typeName);
120  diffusion_ = dict.lookupOrDefault<Switch>("diffusion", false);
121  tolerance_ = dict.lookupOrDefault<scalar>("tolerance", 1e-5);
122 
123  return true;
124 }
125 
126 
128 {
129  return wordList{phiName_, rhoName_};
130 }
131 
132 
134 {
136  (
137  new volScalarField
138  (
139  IOobject
140  (
141  typeName,
142  time_.name(),
143  mesh_,
146  false
147  ),
148  mesh_,
150  patchTypes()
151  )
152  );
153  volScalarField& age = tage.ref();
154 
155  const word divScheme("div(phi," + schemesField_ + ")");
156 
157  // Set under-relaxation coeff
158  scalar relaxCoeff = 0.0;
159  if (mesh_.solution().relaxEquation(schemesField_))
160  {
161  relaxCoeff = mesh_.solution().equationRelaxationFactor(schemesField_);
162  }
163 
166  (
168  );
169 
170  const surfaceScalarField& phi =
171  mesh_.lookupObject<surfaceScalarField>(phiName_);
172 
173  if (phi.dimensions() == dimMass/dimTime)
174  {
175  const volScalarField& rho =
176  mesh_.lookupObject<volScalarField>(rhoName_);
177 
178  const word laplacianScheme("laplacian(muEff," + schemesField_ + ")");
179 
180  tmp<volScalarField> tnuEff;
181  if (diffusion_)
182  {
183  tnuEff = mesh_.lookupType<momentumTransportModel>().nuEff();
184  }
185 
186  for (int i=0; i<=nCorr_; i++)
187  {
188  fvScalarMatrix ageEqn
189  (
190  fvm::div(phi, age, divScheme) == rho + fvModels.source(rho, age)
191  );
192 
193  if (diffusion_)
194  {
195  ageEqn -= fvm::laplacian(rho*tnuEff(), age, laplacianScheme);
196  }
197 
198  ageEqn.relax(relaxCoeff);
199 
200  fvConstraints.constrain(ageEqn);
201 
202  if (converged(i, ageEqn.solve(schemesField_).initialResidual()))
203  {
204  break;
205  };
206 
208  }
209  }
210  else
211  {
212  tmp<volScalarField> tnuEff;
213  word laplacianScheme;
214 
215  if (diffusion_)
216  {
217  tnuEff = mesh_.lookupType<momentumTransportModel>().nuEff();
218 
219  laplacianScheme =
220  "laplacian(" + tnuEff().name() + ',' + schemesField_ + ")";
221  }
222 
223  for (int i=0; i<=nCorr_; i++)
224  {
225  fvScalarMatrix ageEqn
226  (
227  fvm::div(phi, age, divScheme)
229  );
230 
231  if (diffusion_)
232  {
233  ageEqn -= fvm::laplacian(tnuEff(), age, laplacianScheme);
234  }
235 
236  ageEqn.relax(relaxCoeff);
237 
238  fvConstraints.constrain(ageEqn);
239 
240  if (converged(i, ageEqn.solve(schemesField_).initialResidual()))
241  {
242  break;
243  }
244 
246  }
247  }
248 
249  Info<< "Min/max age:" << min(age).value()
250  << ' ' << max(age).value() << endl;
251 
252  return store(tage);
253 }
254 
255 
257 {
258  return writeObject(typeName);
259 }
260 
261 
262 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
const dimensionSet & dimensions() const
Return dimensions.
static const char *const typeName
Definition: Field.H:106
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:100
Abstract base-class for Time/database functionObjects.
Calculates and writes out the time taken for a particle to travel from an inlet to the location....
Definition: age.H:132
age(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: age.C:94
virtual wordList fields() const
Return the list of fields required.
Definition: age.C:127
virtual ~age()
Destructor.
Definition: age.C:108
virtual bool execute()
Execute.
Definition: age.C:133
virtual bool write()
Write.
Definition: age.C:256
virtual bool read(const dictionary &)
Read the data.
Definition: age.C:114
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
const fvMesh & mesh_
Reference to the fvMesh.
Finite volume constraints.
Definition: fvConstraints.H:67
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:603
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:857
Finite volume models.
Definition: fvModels.H:65
tmp< fvMatrix< Type > > source(const VolField< Type > &field) const
Return source for an equation.
Abstract base class for turbulence models (RAS, LES and laminar).
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
Foam::fvConstraints & fvConstraints(Foam::fvConstraints::New(mesh))
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
const scalar nuEff
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for the laplacian of the field.
label patchi
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
Definition: fvmDiv.C:48
Namespace for OpenFOAM.
List< word > wordList
A List of words.
Definition: fileName.H:54
const doubleScalar e
Definition: doubleScalar.H:106
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
messageStream Info
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimTime
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimMass
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
wordList patchTypes(nPatches)
dictionary dict