mixtureFractionSoot.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) 2013-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 
26 #include "mixtureFractionSoot.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 
32 template<class ThermoType>
35 (
36  const fluidThermo& thermo
37 )
38 {
39  if (isA<singleStepReactingMixture<ThermoType> >(thermo))
40  {
41  return dynamic_cast<const singleStepReactingMixture<ThermoType>& >
42  (
43  thermo
44  );
45  }
46  else
47  {
49  (
50  "template<class ThermoType> "
51  "Foam::radiation::mixtureFractionSoot "
52  "("
53  "const dictionary&, "
54  "const fvMesh&"
55  ")"
56  )
57  << "Inconsistent thermo package for " << thermo.type()
58  << "Please select a thermo package based on "
59  << "singleStepReactingMixture" << exit(FatalError);
60 
61  return dynamic_cast<const singleStepReactingMixture<ThermoType>& >
62  (
63  thermo
64  );
65  }
66 
67 }
68 
69 
70 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
71 
72 template<class ThermoType>
74 (
75  const dictionary& dict,
76  const fvMesh& mesh,
77  const word& modelType
78 )
79 :
80  sootModel(dict, mesh, modelType),
81  soot_
82  (
83  IOobject
84  (
85  "soot",
86  mesh_.time().timeName(),
87  mesh_,
88  IOobject::MUST_READ,
89  IOobject::AUTO_WRITE
90  ),
91  mesh_
92  ),
93  coeffsDict_(dict.subOrEmptyDict(modelType + "Coeffs")),
94  nuSoot_(readScalar(coeffsDict_.lookup("nuSoot"))),
95  Wsoot_(readScalar(coeffsDict_.lookup("Wsoot"))),
96  sootMax_(-1),
97  mappingFieldName_
98  (
99  coeffsDict_.lookupOrDefault<word>("mappingFieldName", "none")
100  ),
101  mapFieldMax_(1),
103  mixture_(checkThermo(thermo_))
104 {
105  const Reaction<ThermoType>& reaction = mixture_.operator[](0);
106 
107  const scalarList& specieStoichCoeffs(mixture_.specieStoichCoeffs());
108 
109  scalar totalMol = 0.0;
110  forAll(reaction.rhs(), i)
111  {
112  label speciei = reaction.rhs()[i].index;
113  totalMol += mag(specieStoichCoeffs[speciei]);
114  }
115 
116  totalMol += nuSoot_;
117 
118  scalarList Xi(reaction.rhs().size());
119 
120  scalar Wm = 0.0;
121  forAll(reaction.rhs(), i)
122  {
123  const label speciei = reaction.rhs()[i].index;
124  Xi[i] = mag(specieStoichCoeffs[speciei])/totalMol;
125  Wm += Xi[i]*mixture_.speciesData()[speciei].W();
126  }
127 
128  const scalar XSoot = nuSoot_/totalMol;
129  Wm += XSoot*Wsoot_;
130 
131  sootMax_ = XSoot*Wsoot_/Wm;
132 
133  Info << "Maximum soot mass concentrations: " << sootMax_ << nl;
134 
135  if (mappingFieldName_ == "none")
136  {
137  const label index = reaction.rhs()[0].index;
138  mappingFieldName_ = mixture_.Y(index).name();
139  }
140 
141  const label mapFieldIndex = mixture_.species()[mappingFieldName_];
142 
143  mapFieldMax_ = mixture_.Yprod0()[mapFieldIndex];
144 
145 }
146 
147 
148 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
149 
150 template<class ThermoType>
152 {}
153 
154 
155 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156 
157 template<class ThermoType>
159 {
160  const volScalarField& mapField =
161  mesh_.lookupObject<volScalarField>(mappingFieldName_);
162 
163  soot_ = sootMax_*(mapField/mapFieldMax_);
164 }
165 
166 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
#define readScalar
Definition: doubleScalar.C:38
dimensioned< scalar > mag(const dimensioned< Type > &)
mixtureFractionSoot(const dictionary &dict, const fvMesh &mesh, const word &modelType)
Construct from components.
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
messageStream Info
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const List< specieCoeffs > & rhs() const
Definition: ReactionI.H:59
const word dictName("particleTrackDict")
static const char nl
Definition: Ostream.H:260
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual void correct()
Main update/correction routine.
#define forAll(list, i)
Definition: UList.H:421
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:49
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Simple extension of ReactionThermo to handle reaction kinetics in addition to the equilibrium thermod...
Definition: Reaction.H:53
error FatalError
Info<< "Creating reaction model\n"<< endl;autoPtr< combustionModels::psiCombustionModel > reaction(combustionModels::psiCombustionModel::New(mesh))
Single step reacting mixture.
Base class for soor models.
Definition: sootModel.H:53
psiReactionThermo & thermo
Definition: createFields.H:32
This soot model is purely an state model. The ammount of soot produced is determined by a single step...
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
dictionary subOrEmptyDict(const word &, const bool mustRead=false) const
Find and return a sub-dictionary as a copy, or.
Definition: dictionary.C:675