COxidationDiffusionLimitedRate.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-2018 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 
27 #include "mathematicalConstants.H"
28 
29 using namespace Foam::constant;
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
33 template<class CloudType>
35 (
36  const dictionary& dict,
37  CloudType& owner
38 )
39 :
40  SurfaceReactionModel<CloudType>(dict, owner, typeName),
41  Sb_(readScalar(this->coeffDict().lookup("Sb"))),
42  D_(readScalar(this->coeffDict().lookup("D"))),
43  CsLocalId_(-1),
44  O2GlobalId_(owner.composition().carrierId("O2")),
45  CO2GlobalId_(owner.composition().carrierId("CO2")),
46  WC_(0.0),
47  WO2_(0.0),
48  HcCO2_(0.0)
49 {
50  // Determine Cs ids
51  label idSolid = owner.composition().idSolid();
52  CsLocalId_ = owner.composition().localId(idSolid, "C");
53 
54  // Set local copies of thermo properties
55  WO2_ = owner.thermo().carrier().W(O2GlobalId_);
56  const scalar WCO2 = owner.thermo().carrier().W(CO2GlobalId_);
57  WC_ = WCO2 - WO2_;
58 
59  HcCO2_ = owner.thermo().carrier().Hc(CO2GlobalId_);
60 
61  if (Sb_ < 0)
62  {
64  << "Stoichiometry of reaction, Sb, must be greater than zero" << nl
65  << exit(FatalError);
66  }
67 
68  const scalar YCloc = owner.composition().Y0(idSolid)[CsLocalId_];
69  const scalar YSolidTot = owner.composition().YMixture0()[idSolid];
70  Info<< " C(s): particle mass fraction = " << YCloc*YSolidTot << endl;
71 }
72 
73 
74 template<class CloudType>
76 (
78 )
79 :
81  Sb_(srm.Sb_),
82  D_(srm.D_),
83  CsLocalId_(srm.CsLocalId_),
84  O2GlobalId_(srm.O2GlobalId_),
85  CO2GlobalId_(srm.CO2GlobalId_),
86  WC_(srm.WC_),
87  WO2_(srm.WO2_),
88  HcCO2_(srm.HcCO2_)
89 {}
90 
91 
92 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
93 
94 template<class CloudType>
97 {}
98 
99 
100 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
101 
102 template<class CloudType>
104 (
105  const scalar dt,
106  const label celli,
107  const scalar d,
108  const scalar T,
109  const scalar Tc,
110  const scalar pc,
111  const scalar rhoc,
112  const scalar mass,
113  const scalarField& YGas,
114  const scalarField& YLiquid,
115  const scalarField& YSolid,
116  const scalarField& YMixture,
117  const scalar N,
118  scalarField& dMassGas,
119  scalarField& dMassLiquid,
120  scalarField& dMassSolid,
121  scalarField& dMassSRCarrier
122 ) const
123 {
124  // Fraction of remaining combustible material
125  const label idSolid = CloudType::parcelType::SLD;
126  const scalar fComb = YMixture[idSolid]*YSolid[CsLocalId_];
127 
128  // Surface combustion active combustible fraction is consumed
129  if (fComb < small)
130  {
131  return 0.0;
132  }
133 
134  const SLGThermo& thermo = this->owner().thermo();
135 
136  // Local mass fraction of O2 in the carrier phase
137  const scalar YO2 = thermo.carrier().Y(O2GlobalId_)[celli];
138 
139  // Change in C mass [kg]
140  scalar dmC = 4.0*mathematical::pi*d*D_*YO2*Tc*rhoc/(Sb_*(T + Tc))*dt;
141 
142  // Limit mass transfer by availability of C
143  dmC = min(mass*fComb, dmC);
144 
145  // Change in O2 mass [kg]
146  const scalar dmO2 = dmC/WC_*Sb_*WO2_;
147 
148  // Mass of newly created CO2 [kg]
149  const scalar dmCO2 = dmC + dmO2;
150 
151  // Update local particle C mass
152  dMassSolid[CsLocalId_] += dmC;
153 
154  // Update carrier O2 and CO2 mass
155  dMassSRCarrier[O2GlobalId_] -= dmO2;
156  dMassSRCarrier[CO2GlobalId_] += dmCO2;
157 
158  const scalar HsC = thermo.solids().properties()[CsLocalId_].Hs(T);
159 
160  // carrier sensible enthalpy exchange handled via change in mass
161 
162  // Heat of reaction [J]
163  return dmC*HsC - dmCO2*HcCO2_;
164 }
165 
166 
167 // ************************************************************************* //
Collection of constants.
dictionary dict
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
COxidationDiffusionLimitedRate(const dictionary &dict, CloudType &owner)
Construct from dictionary.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const PtrList< solidProperties > & properties() const
Return the solidProperties properties.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
rhoReactionThermo & thermo
Definition: createFields.H:28
const solidMixtureProperties & solids() const
Return reference to the global (additional) solids.
Definition: SLGThermo.C:134
Diffusion limited rate surface reaction model for coal parcels. Limited to:
virtual scalar calculate(const scalar dt, const label celli, const scalar d, const scalar T, const scalar Tc, const scalar pc, const scalar rhoc, const scalar mass, const scalarField &YGas, const scalarField &YLiquid, const scalarField &YSolid, const scalarField &YMixture, const scalar N, scalarField &dMassGas, scalarField &dMassLiquid, scalarField &dMassSolid, scalarField &dMassSRCarrier) const
Update surface reactions.
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
stressControl lookup("compactNormalStress") >> compactNormalStress
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: SLGThermo.H:62
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
static const char nl
Definition: Ostream.H:265
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const basicSpecieMixture & carrier() const
Return reference to the gaseous components.
Definition: SLGThermo.C:108
messageStream Info
Templated surface reaction model class.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69