The effect of turbulence on mass transfer in solid fuel combustion: RANS model

(1)

ContentslistsavailableatScienceDirect

Combustion and Flame

journalhomepage:www.elsevier.com/locate/combustflame

The effect of turbulence on mass transfer in solid fuel combustion:

RANS model

Ewa Karchniwy

^a^,^b^,^∗

, Nils Erland L. Haugen

^c

, Adam Klimanek

^a

, Øyvind Langørgen

^c

, Sławomir Sładek

^a

aDepartment of Thermal Technology, Silesian University of Technology, Konarskiego 22, Gliwice 44-100, Poland

bDepartment of Energy and Process Engineering, Norwegian University of Science and Technology, Kolbjørn Hejes vei 1B, Trondheim NO-7491, Norway

cSINTEF Energi A.S., Sem Saelands vei 11, Trondheim 7034, Norway

a rt i c l e i nf o

Article history:

Received 20 August 2020 Revised 14 December 2020 Accepted 28 December 2020

Keywords:

Combustion kinetics Combustion rates Char oxidation Diffusion regime

a b s t r a c t

Inthis paper, akinetic-diffusionsurfacecombustion modelis examined. The model ismodifiedsuch thattwoeffectsofturbulenceareincluded:1)enhancementofthemasstransferduetorelativevelocity betweenparticlesand fluid and 2)reduction ofthemass transferdueto turbulence-inducedparticle clustering.Detailsoftheimplementationarediscussedandtheinfluenceofparameterssuchasair-fuel ratio,particlenumberdensity,particlediameter,turbulenceintensityandcharacteristiclengthscalesare studiedtheoretically.Asimplifiednumericalmodelofacombustionchamberiscreatedtoexplorethe effects ofthecombustion model predictions.Finally,the model isincorporatedintosimulationsofan industrial-scaleboilertoinvestigatetheeffectofturbulenceonthenet surfacereactionrate inareal system.Thestudyshowsthatalthoughonaveragethiseffectisratherminor,thereexistregionsinwhich thecarbonconversionrateiseitherdecreasedorincreasedbyturbulence.

ThisisanopenaccessarticleundertheCCBYlicense(http://creativecommons.org/licenses/by/4.0/)

1. Introduction

Modelling of solid fuels combustion and gasification requires taking intoaccount severalimportant processesoccurring during fuel conversion. A solid fuel particleinjected into a hot environ- mentisfirstheatedupanddried.Inthenextstagedevolatilization starts, whichisa complexdecompositionprocess associatedwith the releaseofmultiplegaseous products.During thelast stage of conversion,theremainingcharisconvertedthroughreactionswith the surrounding gas.Inreality, a distinct separationbetweenthe processes can typically not be distinguished, and the drying and devolatilization, aswell asdevolatilizationandcharsurfacereac- tionsoverlap[1,2],inparticularforlargeparticles.Thedevolatiliza- tionismuchfasterthanthecharconversion,especiallyingasifica- tionsystems,whereslowendothermicreactionsareresponsiblefor thecharconversionrate. Manyparametersaffectthedevolatiliza- tion process leadingto differentvolatilecompositions, totalyield andreaction rate.Arangeof modelshavepreviously beendevel- oped, differingconsiderablybytheircomplexityandaccuracy,see

∗Corresponding author at: Department of Thermal Technology, Silesian Univer- sity of Technology, Konarskiego 22, Gliwice 44-100, Poland.

E-mail address: [email protected] (E. Karchniwy).

[2–5] for more detailed information on the process and its modelling.

Thefinal stageoffuelconversion,i.e.thecharconversionpro- cess,isaffectedby:thediffusionofreactantsfromthesurrounding fluid to the particle surface, diffusion within particle pores, het- erogeneousreactionsatexternalandinternalparticlesurfaces(in- cludingreactant gas adsorption anddesorption), evolutionof the char internal structure of pores, ash inhibition and thermal annealing[6]. Severalapproachestocharconversion modelinghave beenproposedintheliterature. Amongthemostcommonlyused isthekinetic-diffusionsurfacereactionratemodel[7,8]according to whichthe overall reaction ratecan be influenced both by the reactionkineticsandthereactantdiffusion.Thismodelusesglobal kineticsand iscomputationally very efficient butit doesnot ex- plicitlyaccountforprocessessuchasevolutionofthecharintrin- sicsurface area andpore diffusion, nor doesit consider changes inparticlediameteranddensity,variationsintheparticlereactiv- ity[9], thermaldeactivation orash inhibition. A much more de- tailedapproachthatincludesalloftheabove-mentionedprocesses istheCarbonBurnoutKinetics(CBK)modelproposedbyHurtetal.

[10]andfurtherextendedtooxidationandgasificationatelevated pressurebyNiksaetal.[11]andLiu&Niksa[12].TheCBKmodel was developed specifically to correctly predict char burnout and isable tocapturea lower reactivityofchars atthefinal stage of

https://doi.org/10.1016/j.combustﬂame.2020.12.040

(2)

conversion.However,thecomputationalexpensemakesthemodel impractical to use in large-scale simulations [6]. More recently, groupsatStanford UniversityandSINTEFhavedevelopedamodel similartotheCBKmodel[13–15].Thismodelhasamoreaccurate descriptionofthesize anddensityevolutionofthechar,together with a detailedintrinsic reaction mechanism. Annealing is,however,notincludedinthismodel.

Solid fuel combustion in industrial-scale facilities most often occurs under turbulent conditions. From the processes involved in char conversionmentioned above, turbulenceprimarily affects the efficiencyof the reactant transport towards the particle surface. This effect of turbulence has been a subject of several re- cent studies [16–19]. Using Direct Numerical Simulations (DNS) and a simplified case in which a passive scalar (reactant) was consumedisothermally, theauthorsof[16,17]showedthat turbulence might have two effects that counteract each other. Krüger et al. [16] demonstrated that the overall conversion rate can be reducediftheturbulent flowpromotes particleclustering.Thisis related to the rapidoxidizer depletion dueto increasedconcen- tration of particles in the clusters. These studies were extended by Haugen et al. [17] who showed that, in addition to particle clustering, turbulence can also increase the rateof heteroge- neousreactionsthroughvelocityfluctuationsthatintensifythere- actant transfer towardsthe particlesurface. Furthermore,Haugen etal.[17] formulatedamodelthatmodifiesthemasstransferco- efficient toaccount forthetwo effectsofturbulenceandverified themodelagainsttheirDNSresults.Theseinvestigationswerefur- therextendedtomorerealistic,non-isothermalconditions[18]and systemsofpolydisperseparticles[19].

In thecurrentstudy,we focusonthe effectsofturbulenceon the mass transfer from the bulk gas to the particle surface. We discuss the model developed by Haugen etal. [17] and apply it to realistic combustion cases by utilizing the Reynolds Averaged Navier-Stokes (RANS)approach. Bothmain effectsare considered:

the enhancement of mass transfer through velocity ﬂuctuations and the mass transfer rate reduction due to turbulence-induced particleclustering. Westudytheparameters affectingtheprocess and show how the two effects of turbulence inﬂuence the char conversioninajetofparticlesandinanindustrial-scaleboiler.

2. Theory

Thereactantconsumptionrateofafuelparticlecanbedeﬁned asthe normalizedquantity relatingtherateofchangeofparticle mass,mp,anditsinitialmass,m_p_,₀,

=− 1 m_p,0

dmp

dt (1)

Inordertoreducecomplexityoftheanalysiswelimit ourdiscus- sion to the context of char burnout. We apply a simple kinetic- diffusionmodel,withapparent ratekinetics.Itshouldbestressed, however,thattheanalysiscanbeeasilyextendedtomoredetailed models.

2.1. Kinetic-diffusionmodel

One ofthemostfrequently usedapproachesinCFDmodelling of solidfuels combustion andgasiﬁcationisto apply thekinetic- diffusionmodel,givenby

=

π

^d²ppox

mp0

1

1/Rdi f+1/Rkin, (2)

wheredpistheparticlediameter,poxisthepartialpressureofox- idizer,R_{di f} isthereactionrateduetodiffusiondeﬁnedas R_{di f}= C

dp

_T₊_T_p

2

³/4

(3)

Fig. 1. The kinetic-diffusion model for T p= T, ρp= 800 kg / m ³, d p= 500 μm , A = 0 . 002 s / m , E = 79 kJ / mol , C = 5 ·10 ⁻¹²sK ⁻³^/⁴.

andR_kin isthekineticreaction rate,whichisoftenwritteninthe Arrheniusform

R_kin=Aexp

− E RTp

. (4)

Inthe above equationsC isa constant, Tp andT arethe temper- aturesofthe particleandofthegassurrounding theparticle,respectively, A isthe pre-exponential factor,E is theactivation energy, andR is the universal gas constant. The kinetic rateR_kin is theapparentrate,thereforetheintrinsicreactivityandporediffu- sionisalreadyaccountedforinparametersAandE.Themodelcan be extendedto account fortheseeffects explicitly,see forexample [20–22].The kinetic-diffusion model,asgiven inEq.(2),was derived withthe assumptionthat the reaction isﬁrst orderwith respecttotheoxidizer ox,seeSmith[20]fordetails.InFig.1the predictionofthekinetic-diffusionmodelisplottedasafunctionof temperatureforagivencondition.Theeffectsofpurechemicalki- netics(R_{di f}=∞)andpurediffusion(R_kin=∞)arealsoshown.As canbe seen,thechemicalreactions areslowatlow temperatures andlimitstheoverallreactionrate(ZoneI).Athightemperatures, thechemicalreactionsarefast,andtheoverallreactionrateislim- ited bythetransport ofoxidizerto theparticlesurface(Zone III).

BetweenzonesIandIII,anintermediatetemperaturerangeexists (ZoneII)inwhichbothchemicalkineticsanddiffusionareimpor- tantindeterminingtheoverallreactionrate.

2.2. Theeffectofmeangas-particlevelocitydifference

TheconstantC entering Eq.(3)incorporatesall theeffects responsible for mass transfer to the particle surface, i.e. the stoi- chiometryofthereaction,diffusiontotheparticlesurfaceandthe effectofconvection.Chenetal.[23]proposed thatfortheithre- action

C_i=s_iM_C M_i

M RT₀⁷^/⁴

p₀

pShD_i_,₀. (5)

Thisisanextension ofaformuladerivedby Baumetal.[7]fora singleoxidationreaction.Here,s_iistheratioofthestoichiometric coeﬃcients ofcarbon and reactant (e.g. s_i=1 for C+O₂→CO₂; s_i=2for2C+O₂→2CO),M_C andM_iarethemolecularweightsof carbonandthereactant ofreaction i,respectively,Misthe mean molecular weight ofthegas inthe particleboundary layer, Shis theSherwoodnumber,D_iisthediffusioncoeﬃcientofthegaseous

(3)

reactantofreactioni,pispressureandsubscript0denotestheref- erencestate.Assumingthattheparticlescanbetreatedasspheres, theSherwoodnumbercanbedeterminedfromtheRanz-Marshall formula[24]

Sh=2.0+0.6Re¹_p^/²Sc¹^/³, (6) where Sc=

ν

/D is the Schmidt number, and Re_p is the particle Reynoldsnumberdeﬁnedas

Rep=

|

^u^p⁻^u

|

^d^p

ν

^, ⁽⁷⁾

where u isthe meangas velocity, up is theparticle velocity and

ν

îs ^the ^kinematic ^viscosity. Ît ^has ^frequently ^been ârgued ^[20]^,

that forfine pulverized fuelparticles the relative particle-gasve- locityissmall,andthusRep→0andSh→2.However,theparticle Reynoldsnumberscanbecomehigherforpressurizedsystemssuch as entrained-flow gasificationreactors [25].Also, forlarger particles,characterizedbylargerStokesnumbers,theeffectcanbecome importantaswell.

Inmodellingofdilute,particulateflows,theeffectofturbulence onparticledispersionisoftenincluded.Oneofthemostfrequently usedapproachesistoapplyastochastictrackingmethod.Insucha casetheparticletrajectoryiscomputedbasedontheinstantaneous fluid velocity, which is a sum ofthe mean fluid velocity andits fluctuatingcomponent,

u=

ζ

2k/3, (8)

where

ζ

îs â ^normally distributed random number and k is the turbulent kinetic energy. Even though this method mayproduce realistic particledispersion,itgivesrisetounphysicallylargerela- tivevelocitydifferencesbetweenparticleandfluid.Rememberthat even tracer particles will experience this unphysical relative velocity, eventhough they in reality willalways followthe fluid in whichtheyareembedded.Thisisbecauseitistheunresolvedtur- bulenteddiesthattransporttheparticles.Inturn,suchanexager- atedrelativevelocitygivestoolargeSherwoodnumberandhence toohightransportrateofmassbetweenfluidandparticle.There- fore, in the following, we use a constant value of the Sherwood numberwhencalculatingCfromEq.(5),andincludetheeffectof turbulencebyapplyingacorrectionfactor

α

˜,aswillbeexplained below.

2.3. Theeffectofturbulenceandparticleclustering

In practical systems the burning particle is exposed to rapid gas velocity ﬂuctuations occurringdue to turbulence. The turbulent motion can be responsible for considerable increase of oxi- dizertransport totheparticlesurfaceduetotheinducedvelocity difference betweentheparticleandthesurroundingﬂuid, asdis- cussedintheformersection.However,theparticlescanalsoform clustersduetoturbulence,whichcanleadtolocaloxidizerdeple- tionandreductionofthereactionrate.Theseeffectswerestudied byDirectNumericalSimulationin[16–19]andthefollowingmodel wasformulatedfortheturbulencecorrectionfactor

α

˜⁼^Sh₂^mod ⁼^Sh₂ _B₊_DaSt^B _/₂

=αcluster

(9)

where ¹₂Shisthe partcorresponding tothe effectof therelative velocity betweenthe particles and the ﬂuid, while

α

cluster is the part that corresponds to clustering. The model parameter B was shownbyHaugenetal.[17]tovarywithStokesnumberas

B=0.08+St/3. (10)

ApartfromtheaverageSherwoodnumber,Sh,two dimensionless numbers enter Eq. (9), namely the Stokes number (St) and the

Damköhlernumber(Da).Theyaredeﬁnedas

St=

τ

^p^/

τ

L, (11)

Da=

τ

L/

τ

c, (12)

where

τ

^p^is^the^particle^response^time,

τ

Listheintegraltimescale ofturbulenceand

τ

^c^is^the^chemical^time^scale,^related^to^the^com-

bustion time.The particleresponse timeis deﬁnedbythe Stokes time

τ

p=

ρ

^p^d²p

18

ρν

^, ⁽¹³⁾

where

ρ

^p^is^the^particle^(material)^density^and

ρ

^is^the^gas^density.

Theintegraltimescalecanbewrittenas

τ

^L⁼²₃^k

ε

⁽¹⁴⁾

andthechemicaltimescaleisdeﬁnedas

1/

τ

^c⁼ⁿ^p^A^p^Sh_D^d^p^, ⁽¹⁵⁾

whereAp=

π

^d²pistheparticleexternalsurfaceareaand np= 6

ρ

^s

ρ

^p

π

^d³p

(16)

istheparticlenumberdensitywith

ρ

^s^being^the^solids^densityⁱⁿ

themixture.Pleasenotethat

ρ

^p^is^the^material^(or^apparent)^den-

sityoftheparticle,whichisvery differentfromthesoliddensity inthemixture,

ρ

s.Forexample,thesoliddensityofcharinairat 1000 K andstoichiometric conditions is around 0.03kg/m³. This corresponds to nearly one hundred 100 μm-sized char particles per cubic centimeter for char particles withan apparent density of600kg/m³.

Inthe following, thephysical reasoningon which Eq.(9)was derived will be described. If the lifetime of a particle cluster (

τ

cluster)isshortrelativetothechemical timescale(

τ

^c^),^the^reac-

tantconcentrationcan beassumedtobe uniformacrosstheclus- terandequaltotheconcentrationoutsidethecluster.Inthiscase, therelevantreactantconsumption rateisgivenbyEq.(15),which is valid for homogeneous distributions of particles and reactant.

For clusters with long lifetimes compared to the chemical time scale,thereactantconcentrationinsidetheclusterisreduced.This means that the overall consumption ratebecomes dependent on clustercharacteristics,suchasclusterdimensionandparticlenum- ber density.The resulting reactant consumption rate(r) is there- forelimitedboth bythe ratedueto Eq.(15)(r_{uni f orm}=1/

τ

^c⁾ ^and

thecluster-characteristic rate(r_cluster=1/

τ

cluster).This means that thereactantconsumptionrate,whichequalsthemasstransferrate, canbewrittenas:

r= r_{uni f orm}r_cluster

r_{uni f orm}+r_cluster. (17)

Therefore, when the above formulation is normalized using the reactant consumption ratefromEq.(15) withSh=2, denotedas r_{uni f orm}

,Sh=2,oneobtains afactorby whichthemasstransferrate isalteredrelativetotheratetypicallyusedinRANSsimulations:

α

˜ ⁼

r_{uni f orm}r_cluster r_{uni f orm}+r_cluster

/r_{uni f orm}_,_Sh₌₂. (18)

Using Eqs.(12) and(15), andaftersome rearranging,Eq.(9) can berecoveredwithB=r_cluster

τ

LSt/2=r_cluster

τ

p/2.Sincer_clusterisun- known,anapproximateexpressionfortheparameterBwasfound usingDNS(see Eq.(10)). Thedetails on theﬁttingprocedure are givenin[17].

Previous studies [16,17] showed that the intensiﬁedtransport ofoxidizertowardstheparticlesurfaceisthedominatingeffectof

(4)

turbulenceatrelativelylowDa.However,astheDamköhlernum- bergetslarger,theimpededreactanttransportassociatedwiththe particleclusteringbecomesthemajorphenomenoncontrollingthe overall surfacereaction rate. It was also found that the effect of clusteringisstrongestwhentheStokes numberisoftheorderof unity.Thereasonforthatisthatsuchconditions(i.e.similarmag- nitudes of particle andﬂow time scales) are themost conducive totheformationofrelativelylong-livedclusters.(Itiswellknown that particle clusters at the Kolmogorov scale, which are due to particleswithKolmogorovbasedStokesnumbersaroundunity,are the strongest and sharpest,but these clusters typically have too shortlifetimestohaveanyrelevanceforthereactanttransport.)

The Sherwoodnumber Sh entering Eq. (9) can still be deter- minedfromEq.(6),however,theparticleReynoldsnumbershould nowbecalculatedas

Rep= ureldp

ν

⁽¹⁹⁾

suchthat theeffectofturbulentvelocityﬂuctuationsistakeninto accountthroughtherelativevelocity,u_rel.Basedonphysicalargu- ments,Haugenetal.[17]proposedthefollowingexpressionforthe averagerelativevelocitydifferencecausedbytheturbulence:

urel=

β

^urms

Stk⁻²_L ^/³−k⁻²_η ^/³

k⁻_L²^/³−k⁻_η²^/³ , (20) where

β

=0.41isamodelconstant,k_Landk_ηaretheintegraland Kolmogorov scale wavenumbers, respectively. The wave numbers canbelinkedtotheturbulentkineticenergyk,itsdissipationrate

ε

,andkinematicviscosity

ν

^as k_L=2

πε

₃

2k

³/2

(21)

k_η=2

π ε ν

³

1/4

. (22)

ThemainassumptionbehindEq.(20)isthattherelativevelocityis induced onlyby thoseturbulent eddiesthat haveturnovertimes,

τ

eddy,that areshorterthan theparticleresponse time,

τ

p. Inthis way,therelativevelocityisproportional tothesquarerootofthe kineticenergy(E(

κ

)⁾^of^thecorrespondingeddies,suchthat:

u_rel∼

_k

eddy

k_η

E

( κ )

^dk

1/2

∼

_k

eddy

k_η

ε

²^/³

κ

⁻⁵^/³^dk

1/2

(23)

where k_ed_d_y=2

π

/

τ

eddyu_ed_d_y. Furthermore,Kolmogorov scaling for the inertial sub-range wasassumed in orderto relate k_ed_d_y with k_L,whilethemodelconstant,

β

,wasobtainedbyfittingthemodel with a large variation of highly accurate direct numerical simulations. It should be mentioned that for very small Stokes numbers thenumeratorof Eq.(20)mightbecome negative. However, at theseconditionsno significantrelative velocity betweenparti- cles andfluid can exist.Therefore, ifthisis thecase, we assume thatu_rel=0.Thiswillresultinatinydiscontinuitiesinthemodel predictionthatwillbevisibleinfigurespresentedinSection3.2.

BycalculatingtheparticleReynoldsnumberbasedontherela- tivevelocityobtainedfromEq.(20),theRanz-Marshallmodel(see Eq.(6))cannowbeusedtoﬁndtheaverageSherwoodnumber,Sh. As canbeseenfromEqs.(19)and(20),theSherwoodnumberSh is affectedbythe turbulenceonly.Thereaction ratedueto diffu- siongivenbyEq.(3)cannowbemodiﬁedtotakeintoaccountthe effectofturbulenceandparticleclusteringas

R_{di f}=

α

˜_d^C

p

_T₊_T

p

2

³/4

. (24)

The model can therefore incorporate the effect of mean gas- particlevelocitythroughEqs.(5)–(7),aswellastheeffectofturbu- lence andparticleclustering throughEqs. (9)–(15).As mentioned

Table 1

The value of γstfor some mixtures.

Mixture Reaction γ^st

Char particles in air C + O 2→ CO 2 11.4 Char particles in 100% CO 2 C + CO 2→ 2 CO 3.7 Char particles in 100% O 2 C + O 2→ CO 2 2.7 Char particles in steam C + H 2O → CO + H 2 1.5 Ilmenite particles in air 4FeTiO 3+ O 2→ 4 TiO 2+ 2 Fe 2O 3 0.225

above,care should betaken whenapplyingEqs. (5)–(7) withthe stochastictrackingmethod.ItshouldalsobestressedthatEqs.(9)–

(15)are suitable to be applied in RANS models,and their form al- lowstodeterminealltherequiredvariablesduringthesimulation.

In thisstudy, themodel wasimplemented into ANSYS Fluent by means ofa User Deﬁned Function(UDF) mechanism. The UDF is providedasasupplementaryﬁletothispaper.

3. Modelsensitivity

Inthissection,numericalexamplesarepresentedinwhichthe modelapplicabilityandtheinfluenceofthemainmodelparame- tersispresented.Thefirst twoexamples arejustgeneralcalcula- tions, whilethe last one is a simplified CFDsimulation. Thisen- ablesustoexaminethepotentialconditionsinwhichtheeffectof turbulencecanbesignificantinpracticalsystems.

3.1. Numericalexample1

Eq.(11) can be re-organized to yieldthe following expression fortheintegraltimescale:

τ

L=

ρ

^p^d²p

18

ρν

^St^. ⁽²⁵⁾

TheDamköhlernumberisthengivenas Da=

τ

L

τ

^c ⁼

τ

L2DnpAp

dp =

ρ

^p^d²p2DnpAp

18

ρν

^St^d^p ^, ⁽²⁶⁾

where,fortheconsiderationsinthissection,ithasbeenassumed thatSh=2.Furthermore,thesolidsdensityinthedomaincanbe expressedas

ρ

^s=

ρ

^g/

γ

^st, (27)

where

ρ

gisthegasdensityinthegas-solidmixtureand

γ

stisthe stoichiometricair-fuelratio.Thevalueof

γ

^st ^for^some ^mixtures^is

giveninTable1.Apartfromchar-basedmixtures,ilmenitewasalso includedduetoitspossibleapplicationsincetheilmeniteparticles canserveasoxygencarriersinChemicalLoopingCombustion.The particlenumberdensitycannowbeexpressedas

np= 6

ρ

^g

π

^d³p

ρ

^p

γ

^st^. ⁽²⁸⁾

Theintrinsicdensityofthegas,

ρ

,is,however,almostthesame asthegaseousdensityofthemixture,

ρ

^g,aslongasthesolidvol- umefractionislow.From theabove,andby usingthat Ap=

π

^d²p, itcanbeshownthat

Da=

ρ

^p^d²p12D

ρ

^g

π

^d²p

18

ρν

^St

π

^d³p

ρ

^p

γ

^st^d^p ⁼

2

ρ

g

3ScSt

γ

st

ρ

^≈³^Sc²^St

γ

st. (29) In Fig. 2, the Damköhler number, as calculated from Eq.(29), is shownasafunctionofStokesnumberforthesamecasesaslisted inTable 1.Clustering isexpectedto slowdown thereactions for Damköhlernumbersaroundorgreaterthanunity[17].FromFig.2, itcanbeseenthatforcarbonoxidationinairtheDamköhlernum- berislargerthan unityonlyforStokesnumberssmallerthan0.1.

(5)

Fig. 2. Damköhler number at stoichiometric conditions as a function of Stokes number for the cases listed in Table 1 . Here, the Schmidt number is set to Sc = 0 . 7 .

Fig. 3. Effect of clustering ( αcluster) as a function of Stokes number for the cases listed in Table 1 .

For oxidation of ilmenite in air, however, the Damköhler number is above 4 even for Stokes number as large as one. Using Eqs.(10)and(29)thepartduetoclusteringcanbeexpressedas

α

cluster= B

B+DaSt/2= 0.08+St/3

0.08+St/3+1/

(

³^Sc

γ

^st

)

^. ⁽³⁰⁾

The value of

α

cluster asa function of Stokes numberis shown in Fig. 3,fromwhich itis clearthat thepotential to reduce the reaction ratehighlydependsonthecomposition ofthemixture.At stoichiometricconditionsthereactionrateduetoclusteringcanbe reducedupto35%inthecaseofcharcombustioninair,whilefor charcombustioninpureO₂orH₂Otheeffectofclusteringcanbe twiceaslarge.Finally,foroxidationofilmeniteinair,thereduction duetoclusteringisdramatic.

Inthisexamplewediscusstheinﬂuenceofselectedmodelpa- rameters on

α

˜^. ^The ^magnitudes ôf ^the ^studied ^parameters ând other essential model parameters are presented in Table 2. They were selected such that they reflect, to some extent, conditions typically found inindustrialscale facilities(reactorsandcombus- tion chambers).Therequired turbulenceparameters,aswouldbe knowninaRANSsimulation,wereestimated.

100⁻⁵ 10⁻⁴ 10⁻³ 10⁻² 0.5

1 1.5 2 2.5 3

˜α,−

I=0.01,L=1m, μ_t/μ=39

100⁻⁵ 10⁻⁴ 10⁻³ 10⁻² 0.5

1 1.5 2 2.5 3

I=0.01,L=10m, μ_t/μ=386

100⁻⁵ 10⁻⁴ 10⁻³ 10⁻² 0.5

1 1.5 2 2.5 3

dp, m

˜α,−

I=0.1,L=1m, μt/μ=386 np=1.0e+04

np=1.0e+05 np=1.0e+06 np=1.0e+07

100⁻⁵ 10⁻⁴ 10⁻³ 10⁻² 0.5

1 1.5 2 2.5 3

dp, m

I=0.1,L=10m, μt/μ=3853 γ= 0.1γ_st

γ=γst γ= 10γst packing limit

Fig. 4. The inﬂuence of parameters from Table 2 on ˜ αfor char particles in air ( ρp = 800 kg / m ³, γst= 11 . 4 ). The legend included in the bottom panels apply to the entire ﬁgure.

Table 2

Studied model input parameters.

Name Symbol Unit Value

Mean gas velocity u m/s 10

Turbulence intensity I – 10 ⁻²; 10 ⁻¹

Domain length scale L m 1; 10

Particle number density n p m ⁻³ 10 ⁴−10 ⁷

Gas density ρ ^kg^/^m³ ^0.35

Gas kinematic viscosity ν ^m²^/^s ¹⁰⁻⁴ Diffusion coeﬃcient D m ²/ s 10 ⁻⁴ Particle (material) density ρp kg / m ³ 800

Inordertocalculatetheturbulencekineticenergykanditsdis- sipation

ε

^the^followingexpressionswereused

urms=uI, (31)

k=3

2u²_rms, (32)

ε

⁼^C³μ^/⁴k³^/²

l , (33)

wherelistheintegrallengthscale,approximatedasl=0.07L[26], andC_μ=0.09[27].Itshouldbenotedthatbyusingsuchadeﬁni- tionofl forlargesystems,theintegrallength scaleislikelytobe overestimated,whichinturnleadstounrealisticallyhighturbulent viscositysince

μ

^t=

ρ

^Cμk²

ε

⁼

ρ

^Cμ¹^/⁴k¹^/²l. (34) Nevertheless,in the absenceof problem-speciﬁcdetails, we stick totheaboveestimation.

Theparticle time scale

τ

^p, thetime scaleof theintegral scale eddies

τ

L and the chemical time scale

τ

c are calculated from Eqs.(13)–(15),respectively.Thesetimescalesarethenusedtocal- culate St and Da, and the mean Sherwoodnumber Sh is calculated using Eqs. (19)–(22) together with Eq. (6). The results, in theformof

α

˜(^d^p)^for^selected^particle^number^densities,âre^presented in Fig. 4 for the case of char particles in air. In the two upperpanels,caseswithlowturbulenceintensity(I=1%)arepre- sented,whereasforthelowerpanelsI=10%.Furthermore,there- sultsshownintheleft-hand sidepanelsdifferfromthose onthe rightside bythe turbulencelength scale,asstatedinthetitle of thefigure.TherearefourblacklinesinFig.4.Thelinewithcircle- shaped markersdividesthe figureintoregions ofrich(below the

(6)

Fig. 5. The inﬂuence of α˜ on the reactant consumption rate for ρp= 800 kg / m ³, d p= 500 μ^m^,^A⁼⁰^.^{002 s}^/^m^,^{E =}^{79 kJ}^/^mol^,^C⁼⁵^·¹⁰⁻¹²^sK⁻³^/⁴^.

line) andlean (abovethe line)conditions, whilethe linewithx- shaped markers corresponds to thepacking limit ofparticles, i.e.

the maximumvolume fractionofthe particles.Forspherical particles, the volume fractionat thepacking limit is assumedto be equalto0.63,whichisatypicallimitforrandomlypacked,spher- icalparticlesofthesamesize.Theremainingtwolinesencompass the region inside which the air-fuel ratio,

γ

, is between0.1 and 10.Itisexpectedthatconditionsinrealsystemscorrespondtothe regionlimitedbythesetwolines.Detailsregardingthederivation ofthepackinglimitlineandthestoichiometriclinecanbefound inAppendixA.

From Fig. 4it is clearthat forthe rangeof examined particle number densities, the turbulence do not have any effect on the mass transfer ifthe particles are too small. This is because particles for which

τ

^p<<

τ

η immediately follow the motion of the ﬂuid,soitisnotpossibleforthemtoformclustersorforthetur- bulence toenhancethe masstransferduetoanyrelativevelocity between ﬂuid and particle. (Pleasenote that for dp of the order of 10⁻⁵−10⁻⁴ Eq.(20)yields negative number inside thesquare root and in this region u_rel=0 was assumed.) For larger particles,whichhavelongerresponsetimes,botheffectsofturbulence canbeobserved.Thelargestmasstransferenhancementis,asex- pected, observed for the high turbulence intensity cases (lower panels of Fig.4), in which

α

˜ ^becomes^greater ^than ¹ ^if^the ^par- ticlenumberdensityissuﬃcientlylow.Forallcasesabove acer- tainnp,theeffectofparticleclusteringbecomesdominant(

α

˜<1).

Thisdecreaseinthereactanttransferrateisparticularlystrongfor the low turbulenceintensitycases(upperpanelsof Fig.4) andit is more intense inlarger facilities (right panels). It isalso worth noticingthat bothscenariosare probablearound thestoichiomet- ric conditions,i.e.wecan expect botheffectsof turbulencetobe observedinrealsystems.

Finally,theinﬂuenceof

α

˜^on^the^reactantconsumptionrate, as givenby Eq. (2) inwhich R_diff is found fromEq.(24), is pre- sentedinFig.5asafunctionoftemperature.Resortingalsotothe resultsshowninFig.4,afactorof2enhancementofreactionrate duetoturbulence(

α

˜>1)canbeexpectedatfavorableﬂowcondi- tionsandhightemperatures.Thereductionoftherate(

α

˜<1)can potentially bemuchstronger.Inthefollowing,wewillinvestigate how

α

˜ ^may^varyⁱⁿ^more^realisticapplications.

In order to visualize and quantifythe effect of turbulenceon pulverizedcharconversion,asimpliﬁedCFDmodelwasdeveloped.

Table 3

Stoichiometric coeﬃcients for reaction (36) and volatiles composition.

C kH lO mN nS o νi

k 1.034 O 2 1.258

l 2.682 CO 1.034

m 0.899 H 2O 1.341 n 0.0274 SO 2 0.0034

o 0.0034 N 2 0.0137

y x

z

jet zero shear stress walls coal + transport air, 10 m/s

8 m

2 m

air co-flow, 5 m/s

Fig. 6. Schematics of the geometry and boundary conditions.

Thegeometryofthemodelwasselectedtobe a2m×2m ×8m cuboid to which coal particles are introduced through a square (4cm × 4cm) inlettogether witha co-ﬂowinghot air. Inside the domainthe particles forma jet andundergodevolatilizationand charcombustion.Themainfeaturesofthenumericalapproachare asfollows.TheNavier-Stokesequationsaresolvedinasteady-state andincompressible form, turbulence ismodelled using the stan- dardk−

ε

^model,^radiationîsâccounted^for^with^the^DiscreteÔr-

dinatesmodelandthe particlesaretrackedinaLagrangian refer- enceframe.Forsimplicity,andsincethefocus ofthepaperison charconversion,thedevolatilizationrateisassumedconstant(=50 1/s).Asinglesurfacereactionisconsidered:

C+O2→CO2 (35)

wherethecorrespondingArrheniusparametersareA=0.002s/m, E=7.9·10⁷J/kmolandthediffusionconstantfromEq.(3)isgiven byC=5·10⁻¹²s/K⁻³^/⁴, while the combustionrate ofvolatiles is computedusingtheFinite-Rate/Eddy-Dissipationmodel,according tothereaction:

C_kH_lOmNnSo+

ν

^O2O2→

ν

^CO^CO⁺

ν

^H2OH2O+

ν

^SO2SO2+

ν

^N2N2, (36) where the stoichiometric coeﬃcients

ν

i and the composition of thefictitious volatilesspeciesC_kH_lO_mN_nS_o are giveninTable 3.A schematicrepresentationofthegeometryandboundaryconditions are given inFig. 6, coal properties are givenin Table 4, andthe mainmodelparametersarepresentedinTable5.Theselectionof thisparticularconfigurationwasmotivatedby thefact that itre- flects typicalconditions forfuel supplyto thecombustion cham-

(7)

Fig. 7. Distribution of Da, St and ˜ αinside the jet.

Table 4 Coal properties.

Proximate analysis Ultimate analysis (daf)

Moisture 0.107 C 0.674

Volatiles 0.446 H 0.05

Fixed carbon 0.357 O 0.267

Ash 0.09 N 0.007

HCV (AR) 22.5 MJ/kg S 0.002

Table 5

CFD model input parameters.

Name Symbol Unit Value

Coal mass flow rate m f kg/s 1 . 5 ·10 ⁻² Transport air mass flow rate m air,1 kg / s 0.0056 Transport air temperature T _air,1 K 1000 Coflow air mass flow rate m _air,2 kg / s 7.0 Coflow air temperature T _air,2 K 1000

Turbulence intensity I – 10 ⁻²

Viscosity ratio μt/μ ^– 50

Coal (material) density ρp kg / m ³ 1400 Coal particle diameter d p m 2 . 5 ·10 ⁻⁴

ber.Moreover,theinputparametersarechosensuchthatthisset- up corresponds (to acertain degree) tothe upper,left-hand side panel ofFig.4,whichmeansthattheturbulenceismostlikelyto reduce themasstransferrate.TherightpanelofFig.7showsthe distributionof

α

˜ ⁱⁿ^a^cross^section^inside^the^jet.^Please^note^that:

1) nointerpolation(nosmoothingbetweencellvalues)isusedto producecontoursof

α

˜ ⁱⁿôrder^toâvoidâ^falseîmpressionôf^low

α

˜ ât^the êdges ôf ^the^jet; ²⁾ ônly^regions ^with^burning^particles aredisplayed.Fromthefigureitcanbeseenthat,fortheconfigu- rationconsidered,theeffectofclusteringissignificant.Infact,

α

˜^is oftheorderof10⁻¹forthemostpartofthejet.Anintensiﬁcation inthemasstransferispredictedonlyattheedgesofthejet,where theparticlenumberdensityislowerandtheturbulenceintensity ishighest. Thereasontheeffectofturbulenceissostrongcanbe understood by inspectingthe DamköhlerandStokes numbers inside the jet. Thesetwo dimensionlessnumbers are showninthe left andmiddlepanelsofFig.7.Eventhough Stdecreasesby 2–3

Fig. 8. Char conversion as a function of distance from the inlet - effect of fuel mass ﬂow rate (in all cases d p= 2 . 5 ·10 ⁻⁴m ).

ordersofmagnitudealongtheparticlejet,Daremains suﬃciently high(of the orderof 1) in the entire volume of the jet to yield

α

˜<1.It shouldalso be notedthat,based onFig. 4,for

α

˜ ^to ^de- creasebelow0.5thelocalconditionsmustcorrespondtoveryrich mixture. Forthe case we study here, a relatively high fuel mass flow rate(m_f=1.5·10⁻²kg/s) waschosen to obtain such conditions butinreality theexistence oflargevolumeswithrich mix- tureisratherunlikelyandmostlyrestrictedtoregionsnexttothe fuelsupply.Therefore,inthefollowingweattempttoverifyifthe effect of turbulence still remains significant for lower fuel mass flowrates.

InFig. 8 charconversion along the jetfor threedifferent fuel massflow rates ispresented. Foreach massflow rate, two cases areshown. Inthefirst,the baselinecase, theeffectofturbulence wasnotaccountedforinthenumericalmodel.Inthesecondcase, theeffectofturbulencewasintroducedthrough theUserDefined Function(UDF).Thiswasdonebymodifying thereactionratedue todiffusionaccordingtoEqs.(9),(10)and(24).Inordertoproduce

(8)

Fig. 9. Contours of temperature (in K) inside the particle jet, from left to right: m f= 1 . 5 ·10 ⁻²kg / s , m f= 1 . 5 ·10 ⁻³kg / s and m f= 1 . 5 ·10 ⁻⁴kg / s .

Fig. 8 the domainwas divided intoN segments along its height.

Foreachsuchsegmentanaveragecarbonconversion(X)wascom- putedforparticleswithinthegivensectionas

X= 1 npart

npart

i=1

Xi= 1 npart

npart

i=1

1− mc

m_c,0

. (37)

In theabove, npart isthe numberof particlespassing through the givensegment, mc andm_c,₀ are thecurrentandinitial parti- clecharmasses,respectively.Theselectedfuelmassﬂowratescan be thoughtofasrich(m_f=1.5·10⁻²kg/s),around-stoichiometric (m_f=1.5·10⁻³kg/s) andlean (m_f=1.5·10⁻⁴kg/s) mixtures, al- though we deliberately do not provide the exact magnitudes of air-fuel ratio(

γ

⁾ âsît ^variessignificantly fromcell to cell. Itcan beseenthatturbulencehasonlyaveryweakpositiveeffectonthe conversionrateifthemassflowrateisverylowor,inotherwords, if

γ

>>

γ

^st^.^The^reason^for^that^is^a^very^low^particle^number^den-

sity, and hence low Damköhler number. In regions with low Da, no dense clusterscan be formed,so thereis noreduction ofthe reactionrateduetoclustering,buta weakincreaseduetoturbu- lence (

α

˜>1)^.^This^behavior îsâlsoⁱⁿâgreement^with^the^results presented inFig.4.Inthecases withhigherfuelmassflow rates inFig.8(redandgreenlines),theparticlenumberdensitiesinthe coreofthejetaremuchhigher,andastrongeffectduetoparticle clusteringcanbeobservedastheconversionismuchslowerinthe caseswherethereactionrateismodifiedbytheUDF.Theconver- sion profiles are similar inboth cases, butthe conversionbegins furtherdownstreamforthecasewiththehighestmassflowrate.

It should be stressed that the results presented in Fig. 8 are strongly affected by the temperature. Even though the same boundary conditions were used, the cases with lower fuel mass flow ratesare characterizedby lower temperaturesinthe system duetothesmalleramountsofreleasedandburnedvolatiles.This isconfirmedinFig.9,wherethecontoursoftemperaturearepre- sented.The consequenceofhighertemperatureishigherreaction rate. This can be observed by comparing the slopes of the conversion profiles in Fig. 8, i.e.the higher the mass flow rate, the

steepertheslope.Atthesametime,asthetemperatureincreases, thediffusionratebecomesmoreimportantintheoverallreaction rate, and thus the observed effect ofturbulence is stronger. The difference in the reaction rates ismore clearly visiblein Fig. 10, which shows contours of the relative rate difference, defined as (−0)/0, where îs^the ^modified^rate încluding^the êffect ofturbulence,and0istheunmodifiedrate.Thehighestrelative ratedifference is observedforthe highestmass flow rate, andit issmaller forthe lower flow rate. It canbe seen that the differ- encesoccurmostlyinthecoreofthejet,wheretheparticlenum- berdensity,andhencetheDamköhlernumber,arethehighest,and thustherateisconsiderablyreduced.However,alsoregionsofin- creasedreactionrateareobservedfurtherawayfromthejet core.

Forthesmallestmass flowratetherelative ratedifferenceis not reducedinthe centerofthejet.Insteada slightlyincreasedreac- tionratecanbeobservedatthejetoutskirts,wheretheturbulence is strongest andthe particle number density is quite low. Based ontheresultsdiscussedabove,wecanconcludethattheeffectof particleclusteringcanbesignificantforaquitewiderangeoffuel mass flowrates or, inother words, fora wide rangeof stoichio- metricconditions.

Another important parameter that inﬂuences how strong the effect of turbulenceis, is the particle size. We examine this pa- rameterby changing theparticle diameter,butkeepingthe mass ﬂow rate constant andequal to 1.5·10⁻³kg/s, corresponding to roughly stoichiometric conditions,forall cases. The particlesizes werechosen suchthat theparticlenumberdensityn_pisdecreas- ingbyafactorof10astheparticlediameterincreases(np∼d⁻³_p ).

AscanbeseeninFig.11,forthesmallestparticles(redlines)the effectofturbulenceamountstoessentiallynodifferencewhenthe total conversion isconsidered (i.e. the distance fromthe inlet to thepointatwhichfullconversionisreached).Ontheotherhand, thelocal reductionofthe conversionrateis actuallyof thesame magnitudeforall particlesizes.Thereasonthedecreasedreaction ratedoesnotaffectthetotalconversiontimeforthesmallestparti- clesisdepletionoftheavailableoxygen,seenasaﬂatteningofthe

(9)

Fig. 10. Contours of relative rate differences (−0)/ 0inside the particle jet, from left to right: m f= 1 . 5 ·10 ⁻²kg / s , m f= 1 . 5 ·10 ⁻³kg / s and m f= 1 . 5 ·10 ⁻⁴kg / s .

Fig. 11. Char conversion as a function of distance from the inlet – effect of particle diameter (in all cases m f= 1 . 5 ·10 ⁻³kg / s ).

’withoutUDF’profileatthefinalstageofconversion.Astheparti- clesizebecomeslarger,particlestravelfurtherdownstreambefore they reachacompleteburnout.Thisisbecausetheselargerparti- clesarelessaffectedbyfluidmotions.Also,onaverage,theyburn in lowertemperatures asthey stillundergo conversionlong after they havepassedthe regions of highesttemperature, i.e.regions of volatile burning. For these larger particles, the effect of turbulence is morepronounced, e.g. particles with dp=2.5·10⁻⁴m (dark bluelines) needto travelaround 50% longertoreach com- pleteburnoutwhentheeffectofturbulenceisaccountedfor.This is opposite to what can be expectedbased on Fig. 4,since fora givenstoichiometriccondition

α

˜ ^increases^for^larger^particle^sizes.

Nevertheless, thedegree to whichthe conversionrateis reduced dependsnotonlyon

α

˜ ^butâlsoôn^the^relative^magnitudesôf^Rkin

andR_{di f}.Astheparticlediameterincreasestheconversionratebe-

comesmorediffusion-controlledsinceR_{di f}∼1/d_p (seeEq.(3)).At thesametime,therateduetokineticsvariesonlyslightly.There- sultingshifttowardsdiffusion-controlledregimeoutweighstheef- fectofhigher

α

˜^and^leads^to^the^conversion^rate^being^reduced^by thesameamount,irrespectiveoftheparticlesize.Finally,itshould be notedthat for evenlarger particles, atsome point

α

˜≥1 (see Fig.4),suchthat noreductionintheconversionrateduetoclus- teringwill be possible, evenfor a fullydiffusion-controlled reaction.Thiswasobserved forparticles withdp∼1·10⁻³mbutwas not shown in Fig.11 dueto much longer time scale required to reachevenafractionalburnout.

Thedegree to whichthe turbulenceinﬂuencesthe surfacereaction rate might also depend on the characteristics of the turbulence itself, such asturbulence intensityor the viscosityratio,

μ

^t/

μ

,astheyare linkedtoturbulencekineticenergyanditsdis- sipation.These twoparameters can affecttheintegral time scale, andthus, the Damköhler andStokes numbers. Their inﬂuence is shownin Fig. 12 from whichit can be seen that the conversion rate is affected in almost exactly the same way for all parameter combinationsthatwe study.Theonlydifference isthatasuf- ﬁciently strong turbulence causes the particles to be converted slightlyfasterasaresultofenhancedmixing.

Finally,weobservedthattheeffectofclusteringweakensifthe jet velocity (velocity atwhich the transport airand particles are introduced) is increased, as shownin Fig. 13. This is due to the Damköhlernumberbeingreducedasthejetvelocityincreases.At evenhigherjetvelocity, theonlyeffectofturbulencewouldbeto increasetheconversionrateasaresultofenhancedmasstransfer totheparticlesurface.

4. Applicationtoanindustrialscaleboiler

In the previous section we explored potential conditions in which turbulence can enhance or decrease the surface reac- tions through the mass transfer rate. As shown by Haugen etal.[17]theseconditionscanbereducedtoonlytwodimension-