wideBandAbsorptionEmission.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-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 \*---------------------------------------------------------------------------*/
25 
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  namespace radiation
34  {
35  defineTypeNameAndDebug(wideBandAbsorptionEmission, 0);
36 
38  (
39  absorptionEmissionModel,
40  wideBandAbsorptionEmission,
41  dictionary
42  );
43  }
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const dictionary& dict,
52  const fvMesh& mesh
53 )
54 :
55  absorptionEmissionModel(dict, mesh),
56  coeffsDict_((dict.optionalSubDict(typeName + "Coeffs"))),
57  speciesNames_(0),
58  specieIndex_(label(0)),
59  lookUpTable_
60  (
61  fileName(coeffsDict_.lookup("lookUpTableFileName")),
62  mesh.time().constant(),
63  mesh
64  ),
66  Yj_(nSpecies_),
67  totalWaveLength_(0)
68 {
69  label nBand = 0;
70  const dictionary& functionDicts = dict.optionalSubDict(typeName +"Coeffs");
71  forAllConstIter(dictionary, functionDicts, iter)
72  {
73  // safety:
74  if (!iter().isDict())
75  {
76  continue;
77  }
78 
79  const dictionary& dict = iter().dict();
80  dict.lookup("bandLimits") >> iBands_[nBand];
81  dict.lookup("EhrrCoeff") >> iEhrrCoeffs_[nBand];
82  totalWaveLength_ += iBands_[nBand][1] - iBands_[nBand][0];
83 
84  label nSpec = 0;
85  const dictionary& specDicts = dict.subDict("species");
86  forAllConstIter(dictionary, specDicts, iter)
87  {
88  const word& key = iter().keyword();
89  if (nBand == 0)
90  {
91  speciesNames_.insert(key, nSpec);
92  }
93  else
94  {
95  if (!speciesNames_.found(key))
96  {
98  << "specie: " << key << "is not in all the bands"
99  << nl << exit(FatalError);
100  }
101  }
102  coeffs_[nBand][nSpec].initialise(specDicts.subDict(key));
103  nSpec++;
104  }
105  nBand++;
106  }
107  nBands_ = nBand;
108 
109  // Check that all the species on the dictionary are present in the
110  // look-up table and save the corresponding indices of the look-up table
111 
112  label j = 0;
113  forAllConstIter(HashTable<label>, speciesNames_, iter)
114  {
115  if (lookUpTable_.found(iter.key()))
116  {
117  label index = lookUpTable_.findFieldIndex(iter.key());
118  Info<< "specie: " << iter.key() << " found in look-up table "
119  << " with index: " << index << endl;
120  specieIndex_[iter()] = index;
121  }
122  else if (mesh.foundObject<volScalarField>(iter.key()))
123  {
124  Yj_.set(j, &mesh.lookupObjectRef<volScalarField>(iter.key()));
125 
126  specieIndex_[iter()] = 0.0;
127  j++;
128  Info<< "species: " << iter.key() << " is being solved" << endl;
129  }
130  else
131  {
133  << "specie: " << iter.key()
134  << " is neither in look-up table : "
135  << lookUpTable_.tableName() << " nor is being solved"
136  << exit(FatalError);
137  }
138  }
139 }
140 
141 
142 
143 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
144 
146 {}
147 
148 
149 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
150 
153 {
154  const volScalarField& T = thermo_.T();
155  const volScalarField& p = thermo_.p();
156  const volScalarField& ft = mesh_.lookupObject<volScalarField>("ft");
157 
158  label nSpecies = speciesNames_.size();
159 
161  (
162  new volScalarField
163  (
164  IOobject
165  (
166  "a",
167  mesh().time().timeName(),
168  mesh(),
171  ),
172  mesh(),
174  )
175  );
176 
177  scalarField& a = ta.ref().primitiveFieldRef();
178 
179  forAll(a, i)
180  {
181  const List<scalar>& species = lookUpTable_.lookUp(ft[i]);
182 
183  for (label n=0; n<nSpecies; n++)
184  {
185  label l = 0;
186  scalar Yipi = 0.0;
187  if (specieIndex_[n] != 0)
188  {
189  // moles x pressure [atm]
190  Yipi = species[specieIndex_[n]]*p[i]*9.869231e-6;
191  }
192  else
193  {
194  // mass fraction from species being solved
195  Yipi = Yj_[l][i];
196  l++;
197  }
198 
199  scalar Ti = T[i];
200 
202  coeffs_[n][bandI].coeffs(T[i]);
203 
204  if (coeffs_[n][bandI].invTemp())
205  {
206  Ti = 1.0/T[i];
207  }
208 
209  a[i]+=
210  Yipi
211  *(
212  ((((b[5]*Ti + b[4])*Ti + b[3])*Ti + b[2])*Ti + b[1])*Ti
213  + b[0]
214  );
215  }
216  }
217 
218  return ta;
219 }
220 
221 
224 {
225  return aCont(bandI);
226 }
227 
228 
231 {
233  (
234  new volScalarField
235  (
236  IOobject
237  (
238  "E",
239  mesh().time().timeName(),
240  mesh(),
243  ),
244  mesh(),
246  )
247  );
248 
249  if (mesh().foundObject<volScalarField>("Qdot"))
250  {
251  const volScalarField& Qdot =
252  mesh().lookupObject<volScalarField>("Qdot");
253 
254  if (Qdot.dimensions() == dimEnergy/dimTime)
255  {
256  E.ref().primitiveFieldRef() =
257  iEhrrCoeffs_[bandI]
258  *Qdot.primitiveField()
259  *(iBands_[bandI][1] - iBands_[bandI][0])
260  /totalWaveLength_
261  /mesh_.V();
262  }
263  else if (Qdot.dimensions() == dimEnergy/dimTime/dimVolume)
264  {
265  E.ref().primitiveFieldRef() =
266  iEhrrCoeffs_[bandI]
267  *Qdot.primitiveField()
268  *(iBands_[bandI][1] - iBands_[bandI][0])
269  /totalWaveLength_;
270  }
271  else
272  {
274  << "Incompatible dimensions for Qdot field" << endl;
275  }
276  }
277 
278  return E;
279 }
280 
281 
283 (
284  volScalarField& a,
285  PtrList<volScalarField>& aLambda
286 ) const
287 {
288  a = dimensionedScalar("zero", dimless/dimLength, 0.0);
289 
290  for (label j=0; j<nBands_; j++)
291  {
292  aLambda[j].primitiveFieldRef() = this->a(j);
293 
294  a.primitiveFieldRef() +=
295  aLambda[j].primitiveField()
296  *(iBands_[j][1] - iBands_[j][0])
297  /totalWaveLength_;
298  }
299 
300 }
301 
302 
303 // ************************************************************************* //
tmp< volScalarField > eCont(const label bandI=0) const
Emission coefficient for continuous phase.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
A class for handling file names.
Definition: fileName.H:69
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Type & lookupObjectRef(const word &name) const
Lookup and return the object reference of the given Type.
defineTypeNameAndDebug(cloudAbsorptionEmission, 0)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
bool foundObject(const word &name) const
Is the named Type found?
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
Model to supply absorption and emission coefficients for radiation modelling.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
scalar Qdot
Definition: solveChemistry.H:2
Macros for easy insertion into run-time selection tables.
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:692
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:759
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
A class for handling words, derived from string.
Definition: word.H:59
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:49
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
void correct(volScalarField &a_, PtrList< volScalarField > &aLambda) const
Correct rays.
word timeName
Definition: getTimeIndex.H:3
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
static const char nl
Definition: Ostream.H:262
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
tmp< volScalarField > ECont(const label bandI=0) const
Emission contribution for continuous phase.
const dimensionSet dimEnergy
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
addToRunTimeSelectionTable(absorptionEmissionModel, cloudAbsorptionEmission, dictionary)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
label n
wideBandAbsorptionEmission(const dictionary &dict, const fvMesh &mesh)
Construct from components.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
tmp< volScalarField > aCont(const label bandI=0) const
Absorption coefficient for continuous phase.
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:576