radiationModel.H
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-2017 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 Namespace
25  Foam::radiation
26 
27 Description
28  Namespace for radiation modelling
29 
30 Class
31  Foam::radiation::radiationModel
32 
33 Description
34  Top level model for radiation modelling
35 
36 SourceFiles
37  radiationModel.C
38  radiationModelNew.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef radiationModel_H
43 #define radiationModel_H
44 
45 #include "IOdictionary.H"
46 #include "autoPtr.H"
47 #include "runTimeSelectionTables.H"
48 #include "volFields.H"
49 #include "fvMatricesFwd.H"
50 #include "Switch.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 
57 class basicThermo;
58 class fvMesh;
59 
60 namespace radiation
61 {
62 
63 // Forward declaration of classes
64 class absorptionEmissionModel;
65 class scatterModel;
66 class sootModel;
67 
68 /*---------------------------------------------------------------------------*\
69  Class radiationModel Declaration
70 \*---------------------------------------------------------------------------*/
71 
72 class radiationModel
73 :
74  public IOdictionary
75 {
76 protected:
77 
78  // Protected data
79 
80  //- Reference to the mesh database
81  const fvMesh& mesh_;
82 
83  //- Reference to the time database
84  const Time& time_;
85 
86  //- Reference to the temperature field
87  const volScalarField& T_;
88 
89  //- Radiation model on/off flag
91 
92  //- Radiation model dictionary
94 
95  //- Radiation solver frequency - number flow solver iterations per
96  // radiation solver iteration
98 
99  //- Flag to enable radiation model to be evaluated on first iteration
100  bool firstIter_;
101 
102 
103  // References to the radiation sub-models
104 
105  //- Absorption/emission model
107 
108  //- Scatter model
110 
111  //- Soot model
113 
114 
115 private:
116 
117  // Private Member Functions
118 
119  //- Create IO object if dictionary is present
120  IOobject createIOobject(const fvMesh& mesh) const;
121 
122  //- Initialise
123  void initialise();
124 
125  //- Disallow default bitwise copy construct
127 
128  //- Disallow default bitwise assignment
129  void operator=(const radiationModel&);
130 
131 
132 public:
133 
134  //- Runtime type information
135  TypeName("radiationModel");
136 
137 
138  // Declare runtime constructor selection table
139 
141  (
142  autoPtr,
144  T,
145  (
146  const volScalarField& T
147  ),
148  (T)
149  );
150 
152  (
153  autoPtr,
155  dictionary,
156  (
157  const dictionary& dict,
158  const volScalarField& T
159  ),
160  (dict, T)
161  );
162 
163 
164  // Constructors
165 
166  //- Null constructor
167  radiationModel(const volScalarField& T);
168 
169  //- Construct from components
170  radiationModel(const word& type, const volScalarField& T);
171 
172  //- Construct from components
174  (
175  const word& type,
176  const dictionary& dict,
177  const volScalarField& T
178  );
179 
180 
181  // Selectors
182 
183  //- Return a reference to the selected radiation model
184  static autoPtr<radiationModel> New(const volScalarField& T);
185 
186  //- Return a reference to the selected radiation model
188  (
189  const dictionary& dict,
190  const volScalarField& T
191  );
192 
193 
194  //- Destructor
195  virtual ~radiationModel();
196 
197 
198  // Member Functions
199 
200  // Edit
201 
202  //- Main update/correction routine
203  virtual void correct();
204 
205  //- Solve radiation equation(s)
206  virtual void calculate() = 0;
207 
208  //- Read radiationProperties dictionary
209  virtual bool read() = 0;
210 
211 
212  // Access
213 
214  //- Radiation model on/off flag
215  const Switch radiation() const
216  {
217  return radiation_;
218  }
219 
220  //- Source term component (for power of T^4)
221  virtual tmp<volScalarField> Rp() const = 0;
222 
223  //- Source term component (constant)
224  virtual tmp<volScalarField::Internal> Ru() const = 0;
225 
226  //- Energy source term
227  virtual tmp<fvScalarMatrix> Sh
228  (
229  const basicThermo& thermo,
230  const volScalarField& he
231  ) const;
232 
233  //- Temperature source term
234  virtual tmp<fvScalarMatrix> ST
235  (
236  const dimensionedScalar& rhoCp,
237  volScalarField& T
238  ) const;
239 
240  //- Access to absorptionEmission model
242 
243  //- Access to soot Model
244  const sootModel& soot() const;
245 };
246 
247 
248 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 #define addToRadiationRunTimeSelectionTables(model) \
251  \
252  addToRunTimeSelectionTable \
253  ( \
254  radiationModel, \
255  model, \
256  dictionary \
257  ); \
258  \
259  addToRunTimeSelectionTable \
260  ( \
261  radiationModel, \
262  model, \
263  T \
264  );
265 
266 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
267 
268 } // End namespace radiation
269 } // End namespace Foam
270 
271 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
272 
273 #endif
274 
275 // ************************************************************************* //
volScalarField & he
Definition: YEEqn.H:51
dictionary dict
dictionary coeffs_
Radiation model dictionary.
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
autoPtr< scatterModel > scatter_
Scatter model.
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:52
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const volScalarField & T_
Reference to the temperature field.
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.
Definition: Switch.H:60
TypeName("radiationModel")
Runtime type information.
const sootModel & soot() const
Access to soot Model.
virtual void calculate()=0
Solve radiation equation(s)
Model to supply absorption and emission coefficients for radiation modelling.
declareRunTimeSelectionTable(autoPtr, radiationModel, T,(const volScalarField &T),(T))
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
static autoPtr< radiationModel > New(const volScalarField &T)
Return a reference to the selected radiation model.
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
bool firstIter_
Flag to enable radiation model to be evaluated on first iteration.
const fvMesh & mesh_
Reference to the mesh database.
const Switch radiation() const
Radiation model on/off flag.
Top level model for radiation modelling.
virtual void correct()
Main update/correction routine.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
virtual tmp< volScalarField::Internal > Ru() const =0
Source term component (constant)
virtual tmp< volScalarField > Rp() const =0
Source term component (for power of T^4)
const absorptionEmissionModel & absorptionEmission() const
Access to absorptionEmission model.
virtual tmp< fvScalarMatrix > Sh(const basicThermo &thermo, const volScalarField &he) const
Energy source term.
Forward declarations of fvMatrix specializations.
const Time & time_
Reference to the time database.
virtual tmp< fvScalarMatrix > ST(const dimensionedScalar &rhoCp, volScalarField &T) const
Temperature source term.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Switch radiation_
Radiation model on/off flag.
Macros to ease declaration of run-time selection tables.
A class for managing temporary objects.
Definition: PtrList.H:53
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:92
autoPtr< absorptionEmissionModel > absorptionEmission_
Absorption/emission model.
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
Namespace for OpenFOAM.