COxidationKineticDiffusionLimitedRate.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-2021 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class CloudType>
34 (
35  const dictionary& dict,
36  CloudType& owner
37 )
38 :
39  SurfaceReactionModel<CloudType>(dict, owner, typeName),
40  Sb_(this->coeffDict().template lookup<scalar>("Sb")),
41  C1_(this->coeffDict().template lookup<scalar>("C1")),
42  C2_(this->coeffDict().template lookup<scalar>("C2")),
43  E_(this->coeffDict().template lookup<scalar>("E")),
44  CsLocalId_(-1),
45  O2GlobalId_(owner.composition().carrierId("O2")),
46  CO2GlobalId_(owner.composition().carrierId("CO2")),
47  WC_(0.0),
48  WO2_(0.0),
49  HcCO2_(0.0)
50 {
51  // Determine Cs ids
52  label idSolid = owner.composition().idSolid();
53  CsLocalId_ = owner.composition().localId(idSolid, "C");
54 
55  // Set local copies of thermo properties
56  WO2_ = owner.composition().carrier().Wi(O2GlobalId_);
57  const scalar WCO2 = owner.composition().carrier().Wi(CO2GlobalId_);
58  WC_ = WCO2 - WO2_;
59 
60  HcCO2_ = owner.composition().carrier().Hf(CO2GlobalId_);
61 
62  const scalar YCloc = owner.composition().Y0(idSolid)[CsLocalId_];
63  const scalar YSolidTot = owner.composition().YMixture0()[idSolid];
64  Info<< " C(s): particle mass fraction = " << YCloc*YSolidTot << endl;
65 }
66 
67 
68 template<class CloudType>
71 (
73 )
74 :
76  Sb_(srm.Sb_),
77  C1_(srm.C1_),
78  C2_(srm.C2_),
79  E_(srm.E_),
80  CsLocalId_(srm.CsLocalId_),
81  O2GlobalId_(srm.O2GlobalId_),
82  CO2GlobalId_(srm.CO2GlobalId_),
83  WC_(srm.WC_),
84  WO2_(srm.WO2_),
85  HcCO2_(srm.HcCO2_)
86 {}
87 
88 
89 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
90 
91 template<class CloudType>
94 {}
95 
96 
97 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
98 
99 template<class CloudType>
101 (
102  const scalar dt,
103  const label celli,
104  const scalar d,
105  const scalar T,
106  const scalar Tc,
107  const scalar pc,
108  const scalar rhoc,
109  const scalar mass,
110  const scalarField& YGas,
111  const scalarField& YLiquid,
112  const scalarField& YSolid,
113  const scalarField& YMixture,
114  const scalar N,
115  scalarField& dMassGas,
116  scalarField& dMassLiquid,
117  scalarField& dMassSolid,
118  scalarField& dMassSRCarrier
119 ) const
120 {
121  // Fraction of remaining combustible material
122  const label idSolid = CloudType::parcelType::SLD;
123  const scalar fComb = YMixture[idSolid]*YSolid[CsLocalId_];
124 
125  // Surface combustion active combustible fraction is consumed
126  if (fComb < small)
127  {
128  return 0.0;
129  }
130 
131  const parcelThermo& thermo = this->owner().thermo();
132  const basicSpecieMixture& carrier = this->owner().composition().carrier();
133 
134  // Local mass fraction of O2 in the carrier phase
135  const scalar YO2 = carrier.Y(O2GlobalId_)[celli];
136 
137  // Diffusion rate coefficient
138  const scalar D0 = C1_/d*pow(0.5*(T + Tc), 0.75);
139 
140  // Kinetic rate
141  const scalar Rk = C2_*exp(-E_/(RR*Tc));
142 
143  // Particle surface area
144  const scalar Ap = constant::mathematical::pi*sqr(d);
145 
146  // Change in C mass [kg]
147  scalar dmC = Ap*rhoc*RR*Tc*YO2/WO2_*D0*Rk/(D0 + Rk)*dt;
148 
149  // Limit mass transfer by availability of C
150  dmC = min(mass*fComb, dmC);
151 
152  // Molar consumption
153  const scalar dOmega = dmC/WC_;
154 
155  // Change in O2 mass [kg]
156  const scalar dmO2 = dOmega*Sb_*WO2_;
157 
158  // Mass of newly created CO2 [kg]
159  const scalar dmCO2 = dOmega*(WC_ + Sb_*WO2_);
160 
161  // Update local particle C mass
162  dMassSolid[CsLocalId_] += dOmega*WC_;
163 
164  // Update carrier O2 and CO2 mass
165  dMassSRCarrier[O2GlobalId_] -= dmO2;
166  dMassSRCarrier[CO2GlobalId_] += dmCO2;
167 
168  const scalar HsC = thermo.solids().properties()[CsLocalId_].Hs(T);
169 
170  // carrier sensible enthalpy exchange handled via change in mass
171 
172  // Heat of reaction [J]
173  return dmC*HsC - dmCO2*HcCO2_;
174 }
175 
176 
177 // ************************************************************************* //
dictionary dict
fluidReactionThermo & thermo
Definition: createFields.H:28
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const PtrList< solidProperties > & properties() const
Return the solidProperties properties.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: parcelThermo.H:58
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Specialisation of basicMixture for a mixture consisting of a number for molecular species...
const dimensionedScalar RR
Universal gas constant: default SI units: [J/kmol/K].
dimensionedScalar exp(const dimensionedScalar &ds)
const solidMixtureProperties & solids() const
Return reference to the global (additional) solids.
Definition: parcelThermo.C:95
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
COxidationKineticDiffusionLimitedRate(const dictionary &dict, CloudType &owner)
Construct from dictionary.
messageStream Info
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
Templated surface reaction model class.
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.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:75
Kinetic/diffusion limited rate surface reaction model for coal parcels. Limited to: ...