greyMeanSolid.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 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 
26 #include "greyMeanSolid.H"
28 #include "unitConversion.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace radiationModels
36 {
37 namespace absorptionEmissionModels
38 {
39  defineTypeNameAndDebug(greyMeanSolid, 0);
40 
42  (
43  absorptionEmissionModel,
44  greyMeanSolid,
45  dictionary
46  );
47 }
48 }
49 }
50 
51 
52 // * * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * //
53 
55 Foam::radiationModels::absorptionEmissionModels::greyMeanSolid::X
56 (
57  const word specie
58 ) const
59 {
60  const volScalarField& T = thermo_.T();
61  const volScalarField& p = thermo_.p();
62 
63  tmp<scalarField> tXj(new scalarField(T.primitiveField().size(), 0.0));
64  scalarField& Xj = tXj.ref();
65 
66  tmp<scalarField> tRhoInv(new scalarField(T.primitiveField().size(), 0.0));
67  scalarField& rhoInv = tRhoInv.ref();
68 
69  forAll(mixture_.Y(), specieI)
70  {
71  const scalarField& Yi = mixture_.Y()[specieI];
72 
73  forAll(rhoInv, iCell)
74  {
75  rhoInv[iCell] +=
76  Yi[iCell]/mixture_.rho(specieI, p[iCell], T[iCell]);
77  }
78  }
79  const scalarField& Yj = mixture_.Y(specie);
80  const label mySpecieI = mixture_.species()[specie];
81  forAll(Xj, iCell)
82  {
83  Xj[iCell] = Yj[iCell]/mixture_.rho(mySpecieI, p[iCell], T[iCell]);
84  }
85 
86  return (Xj/rhoInv);
87 }
88 
89 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
90 
92 (
93  const dictionary& dict,
94  const fvMesh& mesh
95 )
96 :
97  absorptionEmissionModel(dict, mesh),
98  coeffsDict_((dict.optionalSubDict(typeName + "Coeffs"))),
100  speciesNames_(0),
101  mixture_(dynamic_cast<const basicSpecieMixture&>(thermo_)),
102  solidData_(mixture_.Y().size())
103 {
104  if (!isA<basicSpecieMixture>(thermo_))
105  {
107  << "Model requires a multi-component thermo package"
108  << abort(FatalError);
109  }
110 
111  label nFunc = 0;
112  const dictionary& functionDicts = dict.optionalSubDict(typeName + "Coeffs");
113 
114  forAllConstIter(dictionary, functionDicts, iter)
115  {
116  // safety:
117  if (!iter().isDict())
118  {
119  continue;
120  }
121  const word& key = iter().keyword();
122  if (!mixture_.contains(key))
123  {
125  << " specie: " << key << " is not found in the solid mixture"
126  << nl
127  << " specie is the mixture are:" << mixture_.species() << nl
128  << nl << endl;
129  }
130  speciesNames_.insert(key, nFunc);
131  const dictionary& dict = iter().dict();
132  dict.lookup("absorptivity") >> solidData_[nFunc][absorptivity];
133  dict.lookup("emissivity") >> solidData_[nFunc][emissivity];
134 
135  nFunc++;
136  }
137 }
138 
139 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
140 
142 {}
143 
144 
145 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
146 
148 Foam::radiationModels::absorptionEmissionModels::greyMeanSolid::calc
149 (
150  const label propertyId
151 ) const
152 {
154  (
155  new volScalarField
156  (
157  IOobject
158  (
159  "a",
160  mesh().time().timeName(),
161  mesh(),
164  ),
165  mesh(),
167  extrapolatedCalculatedFvPatchVectorField::typeName
168  )
169  );
170 
171  scalarField& a = ta.ref().primitiveFieldRef();
172 
173  forAllConstIter(HashTable<label>, speciesNames_, iter)
174  {
175  if (mixture_.contains(iter.key()))
176  {
177  a += solidData_[iter()][propertyId]*X(iter.key());
178  }
179  }
180 
181  ta.ref().correctBoundaryConditions();
182  return ta;
183 }
184 
185 
188 (
189  const label bandI
190 ) const
191 {
192  return calc(emissivity);
193 }
194 
195 
198 (
199  const label bandI
200 ) const
201 {
202  return calc(absorptivity);
203 }
204 
205 // ************************************************************************* //
addToRunTimeSelectionTable(absorptionEmissionModel, greyMeanCombustion, dictionary)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
tmp< volScalarField > aCont(const label bandI=0) const
Absorption coefficient for continuous phase.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
tmp< volScalarField > eCont(const label bandI=0) const
Emission coefficient for continuous phase.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
Unit conversion functions.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Fundamental solid thermodynamic properties including pressure.
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.
greyMeanSolid(const dictionary &dict, const fvMesh &mesh)
Construct from components.
Definition: greyMeanSolid.C:92
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:123
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
dynamicFvMesh & mesh
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:1001
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
word timeName
Definition: getTimeIndex.H:3
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
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)
static const char nl
Definition: Ostream.H:260
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.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Model to supply absorption and emission coefficients for radiation modelling.
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
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812