radiationModel.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) 2011-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 "radiationModel.H"
28 #include "scatterModel.H"
29 #include "sootModel.H"
30 #include "fvmSup.H"
31 #include "fluidThermo.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  namespace radiation
38  {
39  defineTypeNameAndDebug(radiationModel, 0);
40  defineRunTimeSelectionTable(radiationModel, T);
41  defineRunTimeSelectionTable(radiationModel, dictionary);
42  }
43 }
44 
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47 
48 Foam::IOobject Foam::radiation::radiationModel::createIOobject
49 (
50  const fvMesh& mesh
51 ) const
52 {
53  IOobject io
54  (
55  "radiationProperties",
56  mesh.time().constant(),
57  mesh,
60  );
61 
62  if (io.headerOk())
63  {
65  return io;
66  }
67  else
68  {
69  io.readOpt() = IOobject::NO_READ;
70  return io;
71  }
72 }
73 
74 
75 void Foam::radiation::radiationModel::initialise()
76 {
77  if (radiation_)
78  {
79  solverFreq_ = max(1, lookupOrDefault<label>("solverFreq", 1));
80 
81  absorptionEmission_.reset
82  (
83  absorptionEmissionModel::New(*this, mesh_).ptr()
84  );
85 
86  scatter_.reset(scatterModel::New(*this, mesh_).ptr());
87 
88  soot_.reset(sootModel::New(*this, mesh_).ptr());
89  }
90 }
91 
92 
93 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
94 
95 Foam::radiation::radiationModel::radiationModel(const volScalarField& T)
96 :
98  (
99  IOobject
100  (
101  "radiationProperties",
102  T.time().constant(),
103  T.mesh(),
104  IOobject::NO_READ,
105  IOobject::NO_WRITE
106  )
107  ),
108  mesh_(T.mesh()),
109  time_(T.time()),
110  T_(T),
111  radiation_(false),
112  coeffs_(dictionary::null),
113  solverFreq_(0),
114  firstIter_(true),
115  absorptionEmission_(NULL),
116  scatter_(NULL),
117  soot_(NULL)
118 {}
119 
120 
121 Foam::radiation::radiationModel::radiationModel
122 (
123  const word& type,
124  const volScalarField& T
125 )
126 :
127  IOdictionary(createIOobject(T.mesh())),
128  mesh_(T.mesh()),
129  time_(T.time()),
130  T_(T),
131  radiation_(lookupOrDefault("radiation", true)),
132  coeffs_(subOrEmptyDict(type + "Coeffs")),
133  solverFreq_(1),
134  firstIter_(true),
135  absorptionEmission_(NULL),
136  scatter_(NULL),
137  soot_(NULL)
138 {
139  if (readOpt() == IOobject::NO_READ)
140  {
141  radiation_ = false;
142  }
143 
144  initialise();
145 }
146 
147 
148 Foam::radiation::radiationModel::radiationModel
149 (
150  const word& type,
151  const dictionary& dict,
152  const volScalarField& T
153 )
154 :
156  (
157  IOobject
158  (
159  "radiationProperties",
160  T.time().constant(),
161  T.mesh(),
164  ),
165  dict
166  ),
167  mesh_(T.mesh()),
168  time_(T.time()),
169  T_(T),
170  radiation_(lookupOrDefault("radiation", true)),
171  coeffs_(subOrEmptyDict(type + "Coeffs")),
172  solverFreq_(1),
173  firstIter_(true),
174  absorptionEmission_(NULL),
175  scatter_(NULL),
176  soot_(NULL)
177 {
178  initialise();
179 }
180 
181 
182 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * //
183 
185 {}
186 
187 
188 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
189 
191 {
192  if (regIOobject::read())
193  {
194  lookup("radiation") >> radiation_;
195  coeffs_ = subOrEmptyDict(type() + "Coeffs");
196 
197  solverFreq_ = lookupOrDefault<label>("solverFreq", 1);
198  solverFreq_ = max(1, solverFreq_);
199 
200  return true;
201  }
202  else
203  {
204  return false;
205  }
206 }
207 
208 
210 {
211  if (!radiation_)
212  {
213  return;
214  }
215 
216  if (firstIter_ || (time_.timeIndex() % solverFreq_ == 0))
217  {
218  calculate();
219  firstIter_ = false;
220  }
221 
222  if (!soot_.empty())
223  {
224  soot_->correct();
225  }
226 }
227 
228 
230 (
231  fluidThermo& thermo
232 ) const
233 {
234  volScalarField& he = thermo.he();
235  const volScalarField Cpv(thermo.Cpv());
236  const volScalarField T3(pow3(T_));
237 
238  return
239  (
240  Ru()
241  - fvm::Sp(4.0*Rp()*T3/Cpv, he)
242  - Rp()*T3*(T_ - 4.0*he/Cpv)
243  );
244 }
245 
246 
248 (
249  const dimensionedScalar& rhoCp,
250  volScalarField& T
251 ) const
252 {
253  return
254  (
255  Ru()/rhoCp
256  - fvm::Sp(Rp()*pow3(T)/rhoCp, T)
257  );
258 }
259 
260 
263 {
264  if (!absorptionEmission_.valid())
265  {
267  << "Requested radiation absorptionEmission model, but model is "
268  << "not activate" << abort(FatalError);
269  }
270 
271  return absorptionEmission_();
272 }
273 
274 
277 {
278  if (!soot_.valid())
279  {
281  << "Requested radiation sootModel model, but model is "
282  << "not activate" << abort(FatalError);
283  }
284 
285  return soot_();
286 }
287 
288 
289 // ************************************************************************* //
volScalarField & he
Definition: YEEqn.H:56
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
static autoPtr< scatterModel > New(const dictionary &dict, const fvMesh &mesh)
dictionary coeffs_
Radiation model dictionary.
const sootModel & soot() const
Access to soot Model.
autoPtr< scatterModel > scatter_
Scatter model.
virtual bool read()
Read object.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
dictionary subOrEmptyDict(const word &, const bool mustRead=false) const
Find and return a sub-dictionary as a copy, or.
Definition: dictionary.C:668
const volScalarField & T_
Reference to the temperature field.
static autoPtr< sootModel > New(const dictionary &dict, const fvMesh &mesh)
Selector.
Definition: sootModelNew.C:33
defineTypeNameAndDebug(cloudAbsorptionEmission, 0)
virtual void calculate()=0
Solve radiation equation(s)
const absorptionEmissionModel & absorptionEmission() const
Access to absorptionEmission model.
Model to supply absorption and emission coefficients for radiation modelling.
defineRunTimeSelectionTable(radiationModel, T)
virtual tmp< volScalarField > Cpv() const =0
Heat capacity at constant pressure/volume [J/kg/K].
virtual tmp< fvScalarMatrix > Sh(fluidThermo &thermo) const
Energy source term.
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
virtual ~radiationModel()
Destructor.
dynamicFvMesh & mesh
Base class for soor models.
Definition: sootModel.H:53
virtual bool read()=0
Read radiationProperties dictionary.
autoPtr< sootModel > soot_
Soot model.
A class for handling words, derived from string.
Definition: word.H:59
const Time & time() const
Return time.
Definition: IOobject.C:227
static autoPtr< absorptionEmissionModel > New(const dictionary &dict, const fvMesh &mesh)
Selector.
bool firstIter_
Flag to enable radiation model to be evaluated on first iteration.
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:49
const fvMesh & mesh_
Reference to the mesh database.
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
virtual void correct()
Main update/correction routine.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
virtual tmp< volScalarField > Rp() const =0
Source term component (for power of T^4)
dimensionedScalar pow3(const dimensionedScalar &ds)
virtual tmp< fvScalarMatrix > ST(const dimensionedScalar &rhoCp, volScalarField &T) const
Temperature source term.
IOdictionary(const IOobject &)
Construct given an IOobject.
Definition: IOdictionary.C:45
const Time & time_
Reference to the time database.
virtual tmp< DimensionedField< scalar, volMesh > > Ru() const =0
Source term component (constant)
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
Constant dispersed-phase particle diameter model.
const Mesh & mesh() const
Return mesh.
readOption readOpt() const
Definition: IOobject.H:304
Switch radiation_
Radiation model on/off flag.
A class for managing temporary objects.
Definition: PtrList.H:54
label solverFreq_
Radiation solver frequency - number flow solver iterations per.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
autoPtr< absorptionEmissionModel > absorptionEmission_
Absorption/emission model.
Calculate the matrix for implicit and explicit sources.
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451