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-2019 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 "fvmDdt.H"
28 #include "fvmDiv.H"
30 #include "wallFvPatch.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace functionObjects
39 {
40  defineTypeNameAndDebug(age, 0);
41  addToRunTimeSelectionTable(functionObject, age, dictionary);
42 }
43 }
44 
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47 
48 Foam::wordList Foam::functionObjects::age::patchTypes() const
49 {
50  wordList result
51  (
52  mesh_.boundary().size(),
54  );
55 
56  forAll(mesh_.boundary(), patchi)
57  {
58  if (isA<wallFvPatch>(mesh_.boundary()[patchi]))
59  {
61  }
62  }
63 
64  return result;
65 }
66 
67 
68 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
69 
71 (
72  const word& name,
73  const Time& runTime,
74  const dictionary& dict
75 )
76 :
77  fvMeshFunctionObject(name, runTime, dict),
78  phiName_(),
79  rhoName_(),
80  nCorr_(0),
81  schemesField_()
82 
83 {
84  read(dict);
85 }
86 
87 
88 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
89 
91 {}
92 
93 
94 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
95 
97 {
98  phiName_ = dict.lookupOrDefault<word>("phi", "phi");
99  rhoName_ = dict.lookupOrDefault<word>("rho", "rho");
100 
101  dict.readIfPresent("nCorr", nCorr_);
102 
103  schemesField_ = dict.lookupOrDefault<word>("schemesField", typeName);
104 
105  return true;
106 }
107 
108 
110 {
111  return true;
112 }
113 
114 
116 {
118  (
119  IOobject
120  (
121  typeName,
122  mesh_.time().timeName(),
123  mesh_,
126  ),
127  mesh_,
129  patchTypes()
130  );
131 
132  const word divScheme("div(phi," + schemesField_ + ")");
133 
134  // Set under-relaxation coeff
135  scalar relaxCoeff = 0.0;
136  if (mesh_.relaxEquation(schemesField_))
137  {
138  relaxCoeff = mesh_.equationRelaxationFactor(schemesField_);
139  }
140 
141  // This only works because the null constructed inletValue for an
142  // inletOutletFvPatchField is zero. If we needed any other value we would
143  // have to loop over the inletOutlet patches and explicitly set the
144  // inletValues. We would need to change the interface of inletOutlet in
145  // order to do this.
146 
147  const surfaceScalarField& phi =
148  mesh_.lookupObject<surfaceScalarField>(phiName_);
149 
150  if (phi.dimensions() == dimMass/dimTime)
151  {
152  const volScalarField& rho =
153  mesh_.lookupObject<volScalarField>(rhoName_);
154 
155  for (label i = 0; i <= nCorr_; ++ i)
156  {
157  fvScalarMatrix tEqn
158  (
159  fvm::div(phi, t, divScheme) == rho
160  );
161 
162  tEqn.relax(relaxCoeff);
163 
164  tEqn.solve(schemesField_);
165  }
166  }
167  else
168  {
169  for (label i = 0; i <= nCorr_; ++ i)
170  {
171  fvScalarMatrix tEqn
172  (
173  fvm::div(phi, t, divScheme) == dimensionedScalar(1)
174  );
175 
176  tEqn.relax(relaxCoeff);
177 
178  tEqn.solve(schemesField_);
179  }
180  }
181 
182  Info<< "Min/max age:" << min(t).value() << ' ' << max(t).value() << endl;
183 
184  t.write();
185 
186  return true;
187 }
188 
189 
190 // ************************************************************************* //
virtual bool execute()
Execute.
Definition: age.C:109
virtual bool read(const dictionary &)
Read the data.
Definition: age.C:96
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
surfaceScalarField & phi
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
age(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: age.C:71
addToRunTimeSelectionTable(functionObject, Qdot, dictionary)
wordList patchTypes(nPatches)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
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.
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
Calculate the matrix for the first temporal derivative.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
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)
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:90
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.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
virtual bool write()
Write.
Definition: age.C:115
virtual Ostream & write(const token &)=0
Write next token to stream.
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
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.