FreeFOAM The Cross-Platform CFD Toolkit
greyMeanAbsorptionEmission.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) 2008-2010 OpenCFD Ltd.
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(greyMeanAbsorptionEmission, 0);
36 
38  (
39  absorptionEmissionModel,
40  greyMeanAbsorptionEmission,
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_(0),
59  lookUpTable_
60  (
61  fileName(coeffsDict_.lookup("lookUpTableFileName")),
62  mesh.time().constant(),
63  mesh
64  ),
65  thermo_(mesh.lookupObject<basicThermo>("thermophysicalProperties")),
66  EhrrCoeff_(readScalar(coeffsDict_.lookup("EhrrCoeff"))),
67  Yj_(nSpecies_)
68 {
69  label nFunc = 0;
70  const dictionary& functionDicts = dict.subDict(typeName + "Coeffs");
71 
72  forAllConstIter(dictionary, functionDicts, iter)
73  {
74  // safety:
75  if (!iter().isDict())
76  {
77  continue;
78  }
79  const word& key = iter().keyword();
80  speciesNames_.insert(key, nFunc);
81  const dictionary& dict = iter().dict();
82  coeffs_[nFunc].initialise(dict);
83  nFunc++;
84  }
85 
86  // Check that all the species on the dictionary are present in the
87  // look-up table and save the corresponding indices of the look-up table
88 
89  label j = 0;
90  forAllConstIter(HashTable<label>, speciesNames_, iter)
91  {
92  if (mesh.foundObject<volScalarField>("ft"))
93  {
94  if (lookUpTable_.found(iter.key()))
95  {
96  label index = lookUpTable_.findFieldIndex(iter.key());
97 
98  Info<< "specie: " << iter.key() << " found on look-up table "
99  << " with index: " << index << endl;
100 
101  specieIndex_[iter()] = index;
102  }
103  else if (mesh.foundObject<volScalarField>(iter.key()))
104  {
105  volScalarField& Y =
106  const_cast<volScalarField&>
107  (
108  mesh.lookupObject<volScalarField>(iter.key())
109  );
110  Yj_.set(j, &Y);
111  specieIndex_[iter()] = 0;
112  j++;
113  Info<< "specie: " << iter.key() << " is being solved" << endl;
114  }
115  else
116  {
118  (
119  "Foam::radiation::greyMeanAbsorptionEmission(const"
120  "dictionary& dict, const fvMesh& mesh)"
121  ) << "specie: " << iter.key()
122  << " is neither in look-up table: "
123  << lookUpTable_.tableName()
124  << " nor is being solved" << nl
125  << exit(FatalError);
126  }
127  }
128  else
129  {
131  (
132  "Foam::radiation::greyMeanAbsorptionEmission(const"
133  "dictionary& dict, const fvMesh& mesh)"
134  ) << "specie ft is not present " << nl
135  << exit(FatalError);
136 
137  }
138  }
139 }
140 
141 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
142 
144 {}
145 
146 
147 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
148 
151 {
152  const volScalarField& T = thermo_.T();
153  const volScalarField& p = thermo_.p();
154  const volScalarField& ft = mesh_.lookupObject<volScalarField>("ft");
155 
156  label nSpecies = speciesNames_.size();
157 
159  (
160  new volScalarField
161  (
162  IOobject
163  (
164  "a",
165  mesh().time().timeName(),
166  mesh(),
169  ),
170  mesh(),
172  )
173  );
174 
175  scalarField& a = ta().internalField();
176 
177  forAll(a, i)
178  {
179  const List<scalar>& species = lookUpTable_.lookUp(ft[i]);
180 
181  for (label n=0; n<nSpecies; n++)
182  {
183  label l = 0;
184  scalar Yipi = 0;
185  if (specieIndex_[n] != 0)
186  {
187  //moles x pressure [atm]
188  Yipi = species[specieIndex_[n]]*p[i]*9.869231e-6;
189  }
190  else
191  {
192  // mass fraction
193  Yipi = Yj_[l][i];
194  l++;
195  }
196 
197  const absorptionCoeffs::coeffArray& b = coeffs_[n].coeffs(T[i]);
198 
199  scalar Ti = T[i];
200  // negative temperature exponents
201  if (coeffs_[n].invTemp())
202  {
203  Ti = 1./T[i];
204  }
205  a[i] +=
206  Yipi
207  *(
208  ((((b[5]*Ti + b[4])*Ti + b[3])*Ti + b[2])*Ti + b[1])*Ti
209  + b[0]
210  );
211  }
212  }
213  return ta;
214 }
215 
216 
219 {
221  (
222  new volScalarField
223  (
224  IOobject
225  (
226  "e",
227  mesh().time().timeName(),
228  mesh(),
231  ),
232  mesh(),
234  )
235  );
236 
237  return e;
238 }
239 
240 
243 {
245  (
246  new volScalarField
247  (
248  IOobject
249  (
250  "E",
251  mesh_.time().timeName(),
252  mesh_,
255  ),
256  mesh_,
258  )
259  );
260 
261  if (mesh_.foundObject<volScalarField>("dQ"))
262  {
263  const volScalarField& dQ =
264  mesh_.lookupObject<volScalarField>("dQ");
265  E().internalField() = EhrrCoeff_*dQ;
266  }
267 
268  return E;
269 }
270 
271 
272 // ************************ vim: set sw=4 sts=4 et: ************************ //