CourantNo.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) 2013-2015 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 "CourantNo.H"
27 #include "surfaceFields.H"
28 #include "fvcSurfaceIntegrate.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(CourantNo, 0);
35 }
36 
37 
38 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
39 
41 Foam::CourantNo::byRho
42 (
43  const tmp<volScalarField::DimensionedInternalField>& Co
44 ) const
45 {
46  if (Co().dimensions() == dimDensity)
47  {
48  return Co/obr_.lookupObject<volScalarField>(rhoName_);
49  }
50  else
51  {
52  return Co;
53  }
54 }
55 
56 
57 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
58 
59 Foam::CourantNo::CourantNo
60 (
61  const word& name,
62  const objectRegistry& obr,
63  const dictionary& dict,
64  const bool loadFromFiles
65 )
66 :
67  name_(name),
68  obr_(obr),
69  active_(true),
70  phiName_("phi"),
71  rhoName_("rho")
72 {
73  // Check if the available mesh is an fvMesh, otherwise deactivate
74  if (!isA<fvMesh>(obr_))
75  {
76  active_ = false;
77  WarningIn
78  (
79  "CourantNo::CourantNo"
80  "("
81  "const word&, "
82  "const objectRegistry&, "
83  "const dictionary&, "
84  "const bool"
85  ")"
86  ) << "No fvMesh available, deactivating " << name_ << nl
87  << endl;
88  }
89 
90  read(dict);
91 
92  if (active_)
93  {
94  const fvMesh& mesh = refCast<const fvMesh>(obr_);
95 
96  volScalarField* CourantNoPtr
97  (
98  new volScalarField
99  (
100  IOobject
101  (
102  type(),
103  mesh.time().timeName(),
104  mesh,
107  ),
108  mesh,
109  dimensionedScalar("0", dimless, 0.0),
110  zeroGradientFvPatchScalarField::typeName
111  )
112  );
113 
114  mesh.objectRegistry::store(CourantNoPtr);
115  }
116 }
117 
118 
119 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
120 
122 {}
123 
124 
125 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
126 
128 {
129  if (active_)
130  {
131  phiName_ = dict.lookupOrDefault<word>("phiName", "phi");
132  rhoName_ = dict.lookupOrDefault<word>("rhoName", "rho");
133  }
134 }
135 
136 
138 {
139  if (active_)
140  {
141  const fvMesh& mesh = refCast<const fvMesh>(obr_);
142 
143  const surfaceScalarField& phi =
144  mesh.lookupObject<surfaceScalarField>(phiName_);
145 
146  volScalarField& Co =
147  const_cast<volScalarField&>
148  (
150  );
151 
152  Co.dimensionedInternalField() = byRho
153  (
154  (0.5*mesh.time().deltaT())
155  *fvc::surfaceSum(mag(phi))().dimensionedInternalField()
156  /mesh.V()
157  );
158  Co.correctBoundaryConditions();
159  }
160 }
161 
162 
164 {
165  if (active_)
166  {
167  execute();
168  }
169 }
170 
171 
173 {}
174 
175 
177 {
178  if (active_)
179  {
180  const volScalarField& CourantNo =
181  obr_.lookupObject<volScalarField>(type());
182 
183  Info<< type() << " " << name_ << " output:" << nl
184  << " writing field " << CourantNo.name() << nl
185  << endl;
186 
187  CourantNo.write();
188  }
189 }
190 
191 
192 // ************************************************************************* //
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::surfaceFields.
virtual void read(const dictionary &)
Read the CourantNo data.
Definition: CourantNo.C:127
dimensioned< scalar > mag(const dimensioned< Type > &)
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
tmp< GeometricField< Type, fvPatchField, volMesh > > surfaceSum(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A class for handling words, derived from string.
Definition: word.H:59
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
DimensionedInternalField & dimensionedInternalField()
Return dimensioned internal field.
messageStream Info
dynamicFvMesh & mesh
This function object calculates and outputs the Courant number as a volScalarField. The field is stored on the mesh database so that it can be retrieved and used for other applications.
Definition: CourantNo.H:60
const dimensionSet dimDensity
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
virtual bool write() const
Write using setting from DB.
virtual ~CourantNo()
Destructor.
Definition: CourantNo.C:121
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define WarningIn(functionName)
Report a warning using Foam::Warning.
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField...
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual void write()
Calculate the CourantNo and write.
Definition: CourantNo.C:176
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
virtual void execute()
Execute, currently does nothing.
Definition: CourantNo.C:137
const word & name() const
Return name.
Definition: IOobject.H:260
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
surfaceScalarField & phi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Registry of regIOobjects.
virtual void end()
Execute at the final time-loop, currently does nothing.
Definition: CourantNo.C:163
bool read(const char *, int32_t &)
Definition: int32IO.C:87
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensionedScalar deltaT() const
Return time step.
Definition: TimeState.C:79
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
defineTypeNameAndDebug(combustionModel, 0)
virtual void timeSet()
Called when time was set at the end of the Time::operator++.
Definition: CourantNo.C:172