localDensityAbsorptionEmission.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) 2017 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
31
32// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33
34namespace Foam
35{
36 namespace radiation
37 {
38 defineTypeNameAndDebug(localDensityAbsorptionEmission, 0);
39
41 (
42 absorptionEmissionModel,
43 localDensityAbsorptionEmission,
44 dictionary
45 );
46 }
47}
48
49
50// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51
54{
55 if (!mesh_.foundObject<volScalarField>(alphaName))
56 {
58 << "Unable to retrieve density field " << alphaName << " from "
59 << "database. Available objects:" << mesh_.sortedNames()
60 << exit(FatalError);
61 }
62
63 return mesh_.lookupObject<volScalarField>(alphaName);
64}
65
66
67// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
68
70(
71 const dictionary& dict,
72 const fvMesh& mesh
73)
74:
75 absorptionEmissionModel(dict, mesh),
76 coeffsDict_(dict.subDict(typeName + "Coeffs")),
77 alphaNames_(coeffsDict_.lookup("alphaNames")),
78 aCoeff_(coeffsDict_.lookup("aCoeff")),
79 eCoeff_(coeffsDict_.lookup("eCoeff")),
80 ECoeff_(coeffsDict_.lookup("ECoeff"))
81{}
82
83
84// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
85
88{
89 tmp<volScalarField> ta
90 (
92 (
93 IOobject
94 (
95 "a",
97 mesh_,
100 false
101 ),
102 mesh_,
104 )
105 );
106
107 volScalarField& a = ta.ref();
108
109 forAll(alphaNames_, i)
110 {
111 dimensionedScalar aPhase("a", dimless/dimLength, aCoeff_[i]);
112 a += max(alpha(alphaNames_[i]), scalar(0))*aPhase;
113 }
114
115 return ta;
116}
117
118
121{
122 tmp<volScalarField> te
123 (
125 (
126 IOobject
127 (
128 "e",
129 mesh_.time().timeName(),
130 mesh_,
133 false
134 ),
135 mesh_,
137 )
138 );
139
140 volScalarField& e = te.ref();
141
142 forAll(alphaNames_, i)
143 {
144 dimensionedScalar ePhase("e", dimless/dimLength, eCoeff_[i]);
145 e += max(alpha(alphaNames_[i]), scalar(0))*ePhase;
146 }
147
148 return te;
149}
150
151
154{
155 tmp<volScalarField> tE
156 (
158 (
159 IOobject
160 (
161 "E",
162 mesh_.time().timeName(),
163 mesh_,
166 false
167 ),
168 mesh_,
170 )
171 );
172
173 scalarField& E = tE.ref().primitiveFieldRef();
174
175 forAll(alphaNames_, i)
176 {
177 dimensionedScalar EPhase
178 (
179 "E",
181 ECoeff_[i]
182 );
183
184 E += max(alpha(alphaNames_[i]), scalar(0))*EPhase;
185 }
186
187 return tE;
188}
189
190
191// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
wordList sortedNames() const
The sorted names of all objects.
const Type & lookupObject(const word &name, const bool recursive=false) const
const fvMesh & mesh_
Reference to the fvMesh.
virtual tmp< volScalarField > a(const label bandI=0) const
Absorption coefficient (net)
Constant radiation absorption and emission coefficients for continuous phase.
tmp< volScalarField > eCont(const label bandI=0) const
Emission coefficient for continuous phase.
tmp< volScalarField > ECont(const label bandI=0) const
Emission contribution for continuous phase.
tmp< volScalarField > aCont(const label bandI=0) const
Absorption coefficient for continuous phase.
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
const dimensionSet dimless
Dimensionless.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
volScalarField & alpha
dictionary dict
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333