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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2016-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "radiationModel.H"
31#include "scatterModel.H"
32#include "sootModel.H"
33#include "fvmSup.H"
34#include "basicThermo.H"
35
36// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37
38namespace Foam
39{
40 namespace radiation
41 {
45 }
46}
47
49 "qrExt";
50
52 "qprimaryRad";
53
55 "qreflective";
56
57// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
58
59Foam::IOobject Foam::radiation::radiationModel::createIOobject
60(
61 const fvMesh& mesh
62) const
63{
64 IOobject io
65 (
66 "radiationProperties",
67 mesh.time().constant(),
68 mesh,
71 );
72
73 if (io.typeHeaderOk<IOdictionary>(true))
74 {
76 }
77 else
78 {
80 }
81
82 return io;
83}
84
85
86void Foam::radiation::radiationModel::initialise()
87{
88 if (radiation_)
89 {
90 solverFreq_ = max(1, getOrDefault<label>("solverFreq", 1));
91
92 if (this->found("absorptionEmissionModel"))
93 {
94 absorptionEmission_.reset
95 (
96 absorptionEmissionModel::New(*this, mesh_).ptr()
97 );
98 }
99
100 if (this->found("scatterModel"))
101 {
102 scatter_.reset(scatterModel::New(*this, mesh_).ptr());
103 }
104
105 if (this->found("sootModel"))
106 {
107 soot_.reset(sootModel::New(*this, mesh_).ptr());
108 }
109 }
110}
111
112
113// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
114
116:
118 (
120 (
121 "radiationProperties",
122 T.time().constant(),
123 T.mesh(),
124 IOobject::NO_READ,
125 IOobject::NO_WRITE
126 )
127 ),
128 mesh_(T.mesh()),
129 time_(T.time()),
130 T_(T),
131 radiation_(false),
132 coeffs_(),
133 solverFreq_(0),
134 firstIter_(true),
135 absorptionEmission_(nullptr),
136 scatter_(nullptr),
137 soot_(nullptr)
138{}
139
140
142(
143 const word& type,
144 const volScalarField& T
145)
146:
147 IOdictionary(createIOobject(T.mesh())),
148 mesh_(T.mesh()),
149 time_(T.time()),
150 T_(T),
151 radiation_(getOrDefault("radiation", true)),
152 coeffs_(subOrEmptyDict(type + "Coeffs")),
153 solverFreq_(1),
154 firstIter_(true),
155 absorptionEmission_(nullptr),
156 scatter_(nullptr),
157 soot_(nullptr)
158{
159 if (readOpt() == IOobject::NO_READ)
160 {
161 radiation_ = false;
162 }
163
164 initialise();
165}
166
167
169(
170 const word& type,
171 const dictionary& dict,
172 const volScalarField& T
173)
174:
176 (
178 (
179 "radiationProperties",
180 T.time().constant(),
181 T.mesh(),
182 IOobject::NO_READ,
183 IOobject::NO_WRITE
184 ),
185 dict
186 ),
187 mesh_(T.mesh()),
188 time_(T.time()),
189 T_(T),
190 radiation_(getOrDefault("radiation", true)),
191 coeffs_(subOrEmptyDict(type + "Coeffs")),
192 solverFreq_(1),
193 firstIter_(true),
194 absorptionEmission_(nullptr),
195 scatter_(nullptr),
196 soot_(nullptr)
197{
198 initialise();
199}
200
201
202// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * //
203
205{}
206
207
208// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
209
211{
212 if (regIOobject::read())
213 {
214 readEntry("radiation", radiation_);
215 coeffs_ = subOrEmptyDict(type() + "Coeffs");
216
217 solverFreq_ = getOrDefault<label>("solverFreq", 1);
218 solverFreq_ = max(1, solverFreq_);
219
220 return true;
221 }
222
223 return false;
224}
225
226
228{
229 if (!radiation_)
230 {
231 return;
232 }
233
234 if (firstIter_ || (time_.timeIndex() % solverFreq_ == 0))
235 {
236 calculate();
237 firstIter_ = false;
238 }
239
240 if (soot_)
241 {
242 soot_->correct();
243 }
244}
245
246
248(
249 const basicThermo& thermo,
250 const volScalarField& he
251) const
252{
253 const volScalarField Cpv(thermo.Cpv());
254 const volScalarField T3(pow3(T_));
255
256 return
257 (
258 Ru()
259 - fvm::Sp(4.0*Rp()*T3/Cpv, he)
260 - Rp()*T3*(T_ - 4.0*he/Cpv)
261 );
262}
263
264
266(
269) const
270{
271 return
272 (
273 Ru()/rhoCp
274 - fvm::Sp(Rp()*pow3(T)/rhoCp, T)
275 );
276}
277
278
280(
283) const
284{
285 return
286 (
287 Ru()/rhoCp.ref()
288 - fvm::Sp(Rp()*pow3(T)/rhoCp.ref(), T)
289 );
290}
291
292
294(
296) const
297{
298 return
299 (
300 Ru()
301 - fvm::Sp(Rp()*pow3(T), T)
302 );
303}
304
305
308{
309 if (!absorptionEmission_)
310 {
312 << "Requested radiation absorptionEmission model, but model is "
313 << "not activate" << abort(FatalError);
314 }
315
316 return *absorptionEmission_;
317}
318
319
322{
323 if (!soot_)
324 {
326 << "Requested radiation sootModel model, but model is "
327 << "not activate" << abort(FatalError);
328 }
329
330 return *soot_;
331}
332
333/*
334const Foam::radiation::transmissivityModel&
335Foam::radiation::radiationModel::transmissivity() const
336{
337 if (!transmissivity_)
338 {
339 FatalErrorInFunction
340 << "Requested radiation sootModel model, but model is "
341 << "not activate" << abort(FatalError);
342 }
343
344 return *transmissivity_;
345}
346*/
347
348// ************************************************************************* //
bool found
volScalarField & he
Definition: YEEqn.H:52
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
readOption readOpt() const noexcept
The read option.
Definition: IOobjectI.H:164
bool typeHeaderOk(const bool checkType=true, const bool search=true, const bool verbose=true)
Read header (uses typeFilePath to find file) and check its info.
@ MUST_READ_IF_MODIFIED
Definition: IOobject.H:180
scalar Sh() const
Sherwood number.
const word & constant() const
Return constant name.
Definition: TimePathsI.H:96
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:66
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
constant condensation/saturation model.
Model to supply absorption and emission coefficients for radiation modelling.
Top level model for radiation modelling.
const absorptionEmissionModel & absorptionEmission() const
Access to absorptionEmission model.
virtual ~radiationModel()
Destructor.
static const word externalRadHeatFieldName_
Static name external radiative fluxes.
virtual void correct()
Main update/correction routine.
static const word relfectedFluxName_
Static name for reflected solar fluxes.
virtual bool read()=0
Read radiationProperties dictionary.
Switch radiation_
Radiation model on/off flag.
virtual tmp< fvScalarMatrix > ST(const dimensionedScalar &rhoCp, volScalarField &T) const
Temperature source term.
static const word primaryFluxName_
Static name for primary solar fluxes.
const sootModel & soot() const
Access to soot Model.
Base class for soor models.
Definition: sootModel.H:56
virtual bool read()
Read object.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
const volScalarField & T
dynamicFvMesh & mesh
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Calculate the finiteVolume matrix for implicit and explicit sources.
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
rhoCp
Definition: TEqn.H:3
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensionedScalar pow3(const dimensionedScalar &ds)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
errorManip< error > abort(error &err)
Definition: errorManip.H:144
error FatalError
#define defineRunTimeSelectionTable(baseType, argNames)
Define run-time selection table.
dictionary dict