1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 2008-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*---------------------------------------------------------------------------*/
27 #include "greyMeanAbsorptionEmission.H"
28 #include "addToRunTimeSelectionTable.H"
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 defineTypeNameAndDebug(greyMeanAbsorptionEmission, 0);
38 addToRunTimeSelectionTable
40 absorptionEmissionModel,
41 greyMeanAbsorptionEmission,
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
50 Foam::radiation::greyMeanAbsorptionEmission::greyMeanAbsorptionEmission
52 const dictionary& dict,
56 absorptionEmissionModel(dict, mesh),
57 coeffsDict_((dict.subDict(typeName + "Coeffs"))),
62 fileName(coeffsDict_.lookup("lookUpTableFileName")),
63 mesh.time().constant(),
66 thermo_(mesh.lookupObject<basicThermo>("thermophysicalProperties")),
67 EhrrCoeff_(readScalar(coeffsDict_.lookup("EhrrCoeff"))),
71 const dictionary& functionDicts = dict.subDict(typeName + "Coeffs");
73 forAllConstIter(dictionary, functionDicts, iter)
80 const word& key = iter().keyword();
81 speciesNames_.insert(key, nFunc);
82 const dictionary& dict = iter().dict();
83 coeffs_[nFunc].initialise(dict);
87 // Check that all the species on the dictionary are present in the
88 // look-up table and save the corresponding indices of the look-up table
91 forAllConstIter(HashTable<label>, speciesNames_, iter)
93 if (mesh.foundObject<volScalarField>("ft"))
95 if (lookUpTable_.found(iter.key()))
97 label index = lookUpTable_.findFieldIndex(iter.key());
99 Info<< "specie: " << iter.key() << " found on look-up table "
100 << " with index: " << index << endl;
102 specieIndex_[iter()] = index;
104 else if (mesh.foundObject<volScalarField>(iter.key()))
107 const_cast<volScalarField&>
109 mesh.lookupObject<volScalarField>(iter.key())
112 specieIndex_[iter()] = 0;
114 Info<< "specie: " << iter.key() << " is being solved" << endl;
120 "Foam::radiation::greyMeanAbsorptionEmission(const"
121 "dictionary& dict, const fvMesh& mesh)"
122 ) << "specie: " << iter.key()
123 << " is neither in look-up table: "
124 << lookUpTable_.tableName()
125 << " nor is being solved" << nl
133 "Foam::radiation::greyMeanAbsorptionEmission(const"
134 "dictionary& dict, const fvMesh& mesh)"
135 ) << "specie ft is not present " << nl
142 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
144 Foam::radiation::greyMeanAbsorptionEmission::~greyMeanAbsorptionEmission()
148 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
150 Foam::tmp<Foam::volScalarField>
151 Foam::radiation::greyMeanAbsorptionEmission::aCont(const label bandI) const
153 const volScalarField& T = thermo_.T();
154 const volScalarField& p = thermo_.p();
155 const volScalarField& ft = mesh_.lookupObject<volScalarField>("ft");
157 label nSpecies = speciesNames_.size();
159 tmp<volScalarField> ta
166 mesh().time().timeName(),
172 dimensionedScalar("a", dimless/dimLength, 0.0)
176 scalarField& a = ta().internalField();
180 const List<scalar>& species = lookUpTable_.lookUp(ft[i]);
182 for (label n=0; n<nSpecies; n++)
186 if (specieIndex_[n] != 0)
188 //moles x pressure [atm]
189 Yipi = species[specieIndex_[n]]*p[i]*9.869231e-6;
198 const absorptionCoeffs::coeffArray& b = coeffs_[n].coeffs(T[i]);
201 // negative temperature exponents
202 if (coeffs_[n].invTemp())
209 ((((b[5]*Ti + b[4])*Ti + b[3])*Ti + b[2])*Ti + b[1])*Ti
218 Foam::tmp<Foam::volScalarField>
219 Foam::radiation::greyMeanAbsorptionEmission::eCont(const label bandI) const
221 tmp<volScalarField> e
228 mesh().time().timeName(),
234 dimensionedScalar("e", dimless/dimLength, 0.0)
242 Foam::tmp<Foam::volScalarField>
243 Foam::radiation::greyMeanAbsorptionEmission::ECont(const label bandI) const
245 tmp<volScalarField> E
252 mesh_.time().timeName(),
258 dimensionedScalar("E", dimMass/dimLength/pow3(dimTime), 0.0)
262 if (mesh_.foundObject<volScalarField>("dQ"))
264 const volScalarField& dQ =
265 mesh_.lookupObject<volScalarField>("dQ");
266 E().internalField() = EhrrCoeff_*dQ;
273 // ************************************************************************* //