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-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 
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.subDict(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.subDict(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  (
99  "Foam::radiation::wideBandAbsorptionEmission(const"
100  "dictionary& dict, const fvMesh& mesh)"
101  ) << "specie: " << key << "is not in all the bands"
102  << nl << exit(FatalError);
103  }
104  }
105  coeffs_[nBand][nSpec].initialise(specDicts.subDict(key));
106  nSpec++;
107  }
108  nBand++;
109  }
110  nBands_ = nBand;
111 
112  // Check that all the species on the dictionary are present in the
113  // look-up table and save the corresponding indices of the look-up table
114 
115  label j = 0;
116  forAllConstIter(HashTable<label>, speciesNames_, iter)
117  {
118  if (lookUpTable_.found(iter.key()))
119  {
120  label index = lookUpTable_.findFieldIndex(iter.key());
121  Info<< "specie: " << iter.key() << " found in look-up table "
122  << " with index: " << index << endl;
123  specieIndex_[iter()] = index;
124  }
125  else if (mesh.foundObject<volScalarField>(iter.key()))
126  {
127  volScalarField& Y = const_cast<volScalarField&>
128  (mesh.lookupObject<volScalarField>(iter.key()));
129 
130  Yj_.set(j, &Y);
131 
132  specieIndex_[iter()] = 0.0;
133  j++;
134  Info<< "species: " << iter.key() << " is being solved" << endl;
135  }
136  else
137  {
139  (
140  "radiation::wideBandAbsorptionEmission(const"
141  "dictionary& dict, const fvMesh& mesh)"
142  ) << "specie: " << iter.key()
143  << " is neither in look-up table : "
144  << lookUpTable_.tableName() << " nor is being solved"
145  << exit(FatalError);
146  }
147  }
148 }
149 
150 
151 
152 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
153 
155 {}
156 
157 
158 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
159 
162 {
163  const volScalarField& T = thermo_.T();
164  const volScalarField& p = thermo_.p();
165  const volScalarField& ft = mesh_.lookupObject<volScalarField>("ft");
166 
167  label nSpecies = speciesNames_.size();
168 
170  (
171  new volScalarField
172  (
173  IOobject
174  (
175  "a",
176  mesh().time().timeName(),
177  mesh(),
180  ),
181  mesh(),
183  )
184  );
185 
186  scalarField& a = ta().internalField();
187 
188  forAll(a, i)
189  {
190  const List<scalar>& species = lookUpTable_.lookUp(ft[i]);
191 
192  for (label n=0; n<nSpecies; n++)
193  {
194  label l = 0;
195  scalar Yipi = 0.0;
196  if (specieIndex_[n] != 0)
197  {
198  // moles x pressure [atm]
199  Yipi = species[specieIndex_[n]]*p[i]*9.869231e-6;
200  }
201  else
202  {
203  // mass fraction from species being solved
204  Yipi = Yj_[l][i];
205  l++;
206  }
207 
208  scalar Ti = T[i];
209 
211  coeffs_[n][bandI].coeffs(T[i]);
212 
213  if (coeffs_[n][bandI].invTemp())
214  {
215  Ti = 1.0/T[i];
216  }
217 
218  a[i]+=
219  Yipi
220  *(
221  ((((b[5]*Ti + b[4])*Ti + b[3])*Ti + b[2])*Ti + b[1])*Ti
222  + b[0]
223  );
224  }
225  }
226 
227  return ta;
228 }
229 
230 
233 {
234  return aCont(bandI);
235 }
236 
237 
240 {
242  (
243  new volScalarField
244  (
245  IOobject
246  (
247  "E",
248  mesh().time().timeName(),
249  mesh(),
252  ),
253  mesh(),
255  )
256  );
257 
258  if (mesh().foundObject<volScalarField>("dQ"))
259  {
260  const volScalarField& dQ = mesh().lookupObject<volScalarField>("dQ");
261 
262  if (dQ.dimensions() == dimEnergy/dimTime)
263  {
264  E().internalField() =
265  iEhrrCoeffs_[bandI]
266  *dQ.internalField()
267  *(iBands_[bandI][1] - iBands_[bandI][0])
268  /totalWaveLength_
269  /mesh_.V();
270  }
271  else if (dQ.dimensions() == dimEnergy/dimTime/dimVolume)
272  {
273  E().internalField() =
274  iEhrrCoeffs_[bandI]
275  *dQ.internalField()
276  *(iBands_[bandI][1] - iBands_[bandI][0])
277  /totalWaveLength_;
278  }
279  else
280  {
281  WarningIn
282  (
283  "tmp<volScalarField>"
284  "radiation::wideBandAbsorptionEmission::ECont"
285  "("
286  "const label"
287  ") const"
288  )
289  << "Incompatible dimensions for dQ field" << endl;
290  }
291  }
292 
293  return E;
294 }
295 
296 
298 (
299  volScalarField& a,
300  PtrList<volScalarField>& aLambda
301 ) const
302 {
303  a = dimensionedScalar("zero", dimless/dimLength, 0.0);
304 
305  for (label j=0; j<nBands_; j++)
306  {
307  aLambda[j].internalField() = this->a(j);
308 
309  a.internalField() +=
310  aLambda[j].internalField()
311  *(iBands_[j][1] - iBands_[j][0])
312  /totalWaveLength_;
313  }
314 
315 }
316 
317 
318 // ************************************************************************* //
dimensionedScalar pow3(const dimensionedScalar &ds)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:28
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
bool foundObject(const word &name) const
Is the named Type found?
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
InternalField & internalField()
Return internal field.
const dimensionSet dimEnergy
messageStream Info
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
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
label n
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.
volScalarField & p
Definition: createFields.H:51
PtrList< volScalarField > & Y
Definition: createFields.H:36
#define forAll(list, i)
Definition: UList.H:421
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
dQ
Definition: YEEqn.H:14
Macros for easy insertion into run-time selection tables.
const dimensionSet & dimensions() const
Return dimensions.
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:49
#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
tmp< volScalarField > aCont(const label bandI=0) const
Absorption coefficient for continuous phase.
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
A class for handling file names.
Definition: fileName.H:69
tmp< volScalarField > ECont(const label bandI=0) const
Emission contribution for continuous phase.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
wideBandAbsorptionEmission(const dictionary &dict, const fvMesh &mesh)
Construct from components.
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.
void correct(volScalarField &a_, PtrList< volScalarField > &aLambda) const
Correct rays.
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