greyMeanSolidAbsorptionEmission.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 
28 #include "unitConversion.H"
30 
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  namespace radiation
37  {
38  defineTypeNameAndDebug(greyMeanSolidAbsorptionEmission, 0);
39 
41  (
42  absorptionEmissionModel,
43  greyMeanSolidAbsorptionEmission,
44  dictionary
45  );
46  }
47 }
48 
49 // * * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * //
50 
51 Foam::tmp<Foam::scalarField> Foam::radiation::
52 greyMeanSolidAbsorptionEmission::X(const word specie) const
53 {
54  const volScalarField& T = thermo_.T();
55  const volScalarField& p = thermo_.p();
56 
57  tmp<scalarField> tXj(new scalarField(T.internalField().size(), 0.0));
58  scalarField& Xj = tXj();
59 
60  tmp<scalarField> tRhoInv(new scalarField(T.internalField().size(), 0.0));
61  scalarField& rhoInv = tRhoInv();
62 
63  forAll(mixture_.Y(), specieI)
64  {
65  const scalarField& Yi = mixture_.Y()[specieI];
66 
67  forAll (rhoInv, iCell)
68  {
69  rhoInv[iCell] +=
70  Yi[iCell]/mixture_.rho(specieI, p[iCell], T[iCell]);
71  }
72  }
73  const scalarField& Yj = mixture_.Y(specie);
74  const label mySpecieI = mixture_.species()[specie];
75  forAll(Xj, iCell)
76  {
77  Xj[iCell] = Yj[iCell]/mixture_.rho(mySpecieI, p[iCell], T[iCell]);
78  }
79 
80  return (Xj/rhoInv);
81 }
82 
83 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
84 
87 (
88  const dictionary& dict,
89  const fvMesh& mesh
90 )
91 :
92  absorptionEmissionModel(dict, mesh),
93  coeffsDict_((dict.subDict(typeName + "Coeffs"))),
95  speciesNames_(0),
96  mixture_(dynamic_cast<const basicSpecieMixture&>(thermo_)),
97  solidData_(mixture_.Y().size())
98 {
99  if (!isA<basicSpecieMixture>(thermo_))
100  {
102  (
103  "radiation::greyMeanSolidAbsorptionEmission::"
104  "greyMeanSolidAbsorptionEmission"
105  "("
106  "const dictionary&, "
107  "const fvMesh&"
108  ")"
109  ) << "Model requires a multi-component thermo package"
110  << abort(FatalError);
111  }
112 
113  label nFunc = 0;
114  const dictionary& functionDicts = dict.subDict(typeName + "Coeffs");
115 
116  forAllConstIter(dictionary, functionDicts, iter)
117  {
118  // safety:
119  if (!iter().isDict())
120  {
121  continue;
122  }
123  const word& key = iter().keyword();
124  if (!mixture_.contains(key))
125  {
126  WarningIn
127  (
128  "greyMeanSolidAbsorptionEmission::"
129  "greyMeanSolidAbsorptionEmission "
130  "("
131  " const dictionary& dict,"
132  " const fvMesh& mesh"
133  ")"
134  ) << " specie: " << key << " is not found in the solid mixture"
135  << nl
136  << " specie is the mixture are:" << mixture_.species() << nl
137  << nl << endl;
138  }
139  speciesNames_.insert(key, nFunc);
140  const dictionary& dict = iter().dict();
141  dict.lookup("absorptivity") >> solidData_[nFunc][absorptivity];
142  dict.lookup("emissivity") >> solidData_[nFunc][emissivity];
143 
144  nFunc++;
145  }
146 }
147 
148 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
149 
152 {}
153 
154 
155 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156 
158 Foam::radiation::greyMeanSolidAbsorptionEmission::
159 calc(const label propertyId) const
160 {
162  (
163  new volScalarField
164  (
165  IOobject
166  (
167  "a",
168  mesh().time().timeName(),
169  mesh(),
172  ),
173  mesh(),
175  zeroGradientFvPatchVectorField::typeName
176  )
177  );
178 
179  scalarField& a = ta().internalField();
180 
181  forAllConstIter(HashTable<label>, speciesNames_, iter)
182  {
183  if (mixture_.contains(iter.key()))
184  {
185  a += solidData_[iter()][propertyId]*X(iter.key());
186  }
187  }
188 
189  ta().correctBoundaryConditions();
190  return ta;
191 }
192 
193 
196 (
197  const label bandI
198 ) const
199 {
200  return calc(emissivity);
201 }
202 
203 
206 (
207  const label bandI
208 ) const
209 {
210  return calc(absorptivity);
211 }
212 
213 // ************************************************************************* //
void calc(const argList &args, const Time &runTime, const fvMesh &mesh)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
tmp< volScalarField > aCont(const label bandI=0) const
Absorption coefficient for continuous phase.
defineTypeNameAndDebug(cloudAbsorptionEmission, 0)
A class for handling words, derived from string.
Definition: word.H:59
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
addToRunTimeSelectionTable(absorptionEmissionModel, cloudAbsorptionEmission, dictionary)
dynamicFvMesh & mesh
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:638
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
Namespace for OpenFOAM.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
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.
Model to supply absorption and emission coefficients for radiation modelling.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
#define forAll(list, i)
Definition: UList.H:421
Fundamental solid thermodynamic properties.
Definition: solidThermo.H:54
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Macros for easy insertion into run-time selection tables.
Unit conversion functions.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
greyMeanSolidAbsorptionEmission(const dictionary &dict, const fvMesh &mesh)
Construct from components.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
tmp< volScalarField > eCont(const label bandI=0) const
Emission coefficient for continuous phase.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
word timeName
Definition: getTimeIndex.H:3