Shock Mach number influence on reaction wave types and mixing in reactive shock–bubble interaction

(1)

Shock Mach number influence on reaction wave types and mixing in reactive

shock–bubble interaction

Diegelmann, F; Hickel, S; Adams, NA

DOI

10.1016/j.combustflame.2016.09.014

Publication date

2016

Document Version

Final published version

Published in

Combustion and Flame

Citation (APA)

Diegelmann, F., Hickel, S., & Adams, NA. (2016). Shock Mach number influence on reaction wave types

and mixing in reactive shock–bubble interaction. Combustion and Flame, 174, 85-99.

https://doi.org/10.1016/j.combustflame.2016.09.014

Important note

To cite this publication, please use the final published version (if applicable).

Please check the document version above.

Copyright

Other than for strictly personal use, it is not permitted to download, forward or distribute the text or part of it, without the consent of the author(s) and/or copyright holder(s), unless the work is under an open content license such as Creative Commons. Takedown policy

Please contact us and provide details if you believe this document breaches copyrights. We will remove access to the work immediately and investigate your claim.

This work is downloaded from Delft University of Technology.

(2)

ContentslistsavailableatScienceDirect

Combustion

and

Flame

journalhomepage:www.elsevier.com/locate/combustflame

Shock

Mach

number

inﬂuence

on

reaction

wave

types

and

mixing

in

reactive

shock–bubble

interaction

Felix

Diegelmann

a,∗

_,

_Stefan

_Hickel

a,b

_,

_Nikolaus

_A.

_Adams

a

a Institute of Aerodynamics and Fluid Mechanics, Technische Universität München, Garching 85748, Germany b Faculty of Aerospace Engineering, TU Delft, 2629 HS Delft, Netherlands

a

r

t

i

c

l

e

i

n

f

o

Article history:

Received 4 March 2016 Revised 3 August 2016 Accepted 13 September 2016 Available online 30 September 2016 Keywords: Shock wave Richtmyer–Meshkov instability Shock-bubble interaction Detonation Deﬂagration

a

b

s

t

r

a

c

t

Wepresentnumericalsimulationsforareactiveshock–bubbleinteractionwithdetailedchemistry.The convexshapeofthebubbleleadstoshockfocusing,whichgeneratesspotsofhighpressureand tempera-ture.PressureandtemperaturelevelsaresufficienttoignitethestoichiometricH2–O2gasmixture.Shock MachnumbersbetweenMa=2.13andMa=2.90inducedifferentreactionwavetypes(deflagrationand detonation).Depending onthe shockMach numberlow-pressurereactions orhigh-pressurechemistry areprevalent.AdeflagrationwaveisobservedforthelowestshockMachnumber.ShockMachnumbers ofMa=2.30orhigherignitethegasmixtureafterashortinductiontime,followedbyadetonationwave. AnintermediateshockstrengthofMa=2.19inducesdeflagrationthattransitionsintoadetonationwave. Richtmyer–MeshkovandKelvin–Helmholtzinstabilityevolutionsexhibitahighsensitivitytothereaction wavetype,whichinturnhasdistincteffectsonthespatialandtemporalevolutionofthegasbubble.We observeasignificantreductioninmixingforbothreactionwavetypes,wherein detonationshowsthe strongesteffect.Furthermore,weobserveaverygoodagreementwithexperimentalobservations.

1. Introduction

The interaction between high-speed reactive flows and shock waves is a generic situation present in many combustion sys-tems. Controlled application can promote mixing; uncontrolled interactions, however, can lead to undesirable heat release and thermomechanical loads. Especially in supersonic combustion, where the rapid and efficient mixing of fuel and oxidizer is crucial, astheresidencetime of thefuel–oxidizer mixture inthe combustionchamberisonlyafewmilliseconds [1],mixingcanbe enhanced sufficiently by shock-induced instabilities. The selected genericconfigurationofreactingshock–bubbleinteraction(RSBI)is representative foralarge rangeofhydrodynamicinstabilities and differentreaction wavetypesoccurringinapplication, andallows ustostudytheinteractionbetweendifferenteffectsindetail.

1.1. Hydrodynamicinstabilities

Two hydrodynamic instabilities dominate in a RSBI: the Richtmyer–Meshkov instability (RMI) and the Kelvin–Helmholtz instability (KHI). RMI can enhance mixingin high-speed reactive

∗ _{Corresponding author.}

E-mail addresses: felix.diegelmann@aer.mw.tum.de , felix.diegelmann@gmail.com (F. Diegelmann).

flows, promote turbulent mixing and thus increase the burning efficiencyofsupersoniccombustionengines [2].Theinstability oc-cursattheinterfacebetweentwofluidsofdifferentdensities. The-oreticallystatedin1960byRichtmyer [3]andexperimentally veri-fiedbyMeshkov [4]in1969,RMIcanbeconsideredasthe impul-sivelimitoftheRayleigh–Taylorinstability [5,6].Themisalignment ofpressuregradient,

∇

p,associatedwithashockwaveanddensity gradient,

∇ρ

, atthe material interface causesbaroclinic vorticity productionattheinterface.Forcomprehensivereviewsthereader isreferredtoBrouillette [7]andZabusky [8].RMIoccursonawide rangeofhighlyreactiveenvironmentsfromextremelylargescales inastrophysics [9],tointermediatescalesincombustion [1,10]and downtoverysmallscalesininertialconﬁnementfusion [11].

RMI induces velocity shearandsmall perturbations atthe in-terfaceof the bubble,which are necessary preconditions forKHI

[12]. The perturbations are ampliﬁed, eventually generating vor-tices at the interface accompanied by the appearance of smaller scales [7].KHIdrivesthebreakupoflarge-scalestructures [13]and forcesmixing [14].Botheffectsarethemainhydrodynamicdrivers inRSBI.

1.2.Shock-inducedignitionandreactionwaves

Independently of the scale, RMI is accompanied by a second phenomenoninreactivegasmixtures:theshock-inducedvariation

http://dx.doi.org/10.1016/j.combustﬂame.2016.09.014

0010-2180/© 2016 The Authors. Published by Elsevier Inc. on behalf of The Combustion Institute. This is an open access article under the CC BY-NC-ND license ( http://creativecommons.org/licenses/by-nc-nd/4.0/ ).

(3)

ofthermodynamicproperties,whichcanleadtoignition,followed byareactionwave.Tworeactionwavetypescanbedistinguished: deﬂagration and detonation.Deﬂagration is a subsonic diffusion-drivenreactionwavethatpropagatesthroughthegasmixturedue todirecttransferofchemicalenergyfromburningtounburnedgas

[15].Detonation isdriven bya fastchemical reactionandthe as-sociatedlargeheatreleasewithinthereactionwave.Ashockwave immediatelyprecedes the detonation wave and preheats thegas mixturebycompression [15].The detonationwavepropagatesup to 108 _times _faster _than _the _{deﬂagration} _wave_[16]_. _Due _to _the

largedifferencesinthecharacteristicreactiontime scales,the re-actionwavetypehasacrucialinﬂuenceontheﬂowevolution.

Under certain circumstances a deflagration wave can trans-formintoadetonationwave.Deflagration-to-detonationtransition (DDT)isoneofthemostinterestingunresolvedproblemsin com-bustion theory.Generally, a self-propagating deflagration wave is unstable and tends to accelerate. Under specific conditions the continuousacceleration cansuddenly transitionintoadetonation wave [17]. Liberman et al. [18] proposed a mechanism mainly driven by flame acceleration divided into three stages.The reac-tionfront acceleratesandproducesshockwavesfar aheadofthe flame.Thereafter,theaccelerationdecreases,shocksareformedon theflamesurfaceandpocketsofcompressedandheatedunburnt gasemerge(preheatzone).Inthefinalstagethetransitionto det-onationhappens:theflamepropagatesintothepreheat zoneand producesalargeamplitudepressurepulse.Increasingpressure en-hancesreactionratesandthefeedbackbetweenthepressurepeak andthe reaction leads to a growth of the pressure peak, which steepensintoastrongshockthat,coupledwiththereactionzone, finallyformsanoverdrivendetonationwave.

Furthermore,theflamefrontcanpropagateintoregions ofgas that already have been compressed and preheated by preceding shockwaves such as in shock–bubble interactions (SBI). The re-action ratesand the heat release are enhanced in theseregions, which in turn increases the pressure pulse and accelerates the transitiontodetonation.Ingeneral,DDTcanoccurintworegions: itdevelopsfromthepreheated,compressed gasmixturebetween theleadingshockwaveandtheflameoritarisesfromwithin the flame [19]. The latter transition process is relevant for the pre-sentedstudyasRSBIcontainsregions ofirregularcompressionby theinitialshockwave.

1.3.Reactingshock–bubbleinteraction

Theimpactofashockwaveonareactivegasbubbleallowsto investigate the interaction between shock-induced hydrodynamic instabilities and ignition. The shock wave triggers RMI and the pressureandtemperatureincrease leadstothe formationof rad-icals,whichaccumulateuntilthegasmixtureignites.RMI,dueto themisalignmentofthepressureanddensitygradientatthe bub-bleinterface,causes thebubbleto evolve intoavortex ring. Pro-videdthattheinitialkinetic-energyinputissufficient,theflow de-velopsaturbulentmixingzonethroughnon-linearinteractions of thematerial interface perturbations [7,8].Upon contact, the inci-dentshockwaveispartiallyreflectedandpartiallytransmitted.For anAtwoodnumberA=

(

ρ

1−

ρ

2

)

/

(

ρ

1+

ρ

2

)

<0(thebubblegasis

lighterthanthe ambientgas),thetransmittedshockwave propa-gates fasterthanthe incident shockwave.A _> 0showsthe con-verse effect,the transmitted shockwave travelsslower than the incidentshockwave outsideofthebubble.Thetransmittedshock wavefocusesatthedownstreampoleofthebubbleandcollapses intoasinglepoint(shock-focusingpoint).

Classical inertSBI wasthe subjectof severalstudies overthe lastdecades.HaasandSturtevant [20] investigatedtheinteraction ofshock waves propagating in airwith a gas bubble ﬁlled with eitherheliumorR-22.Theirexperimentalresultscontributedtoa

betterunderstandingofthetemporalbubbleevolutionundershock accelerationandestablished anewclassofcanonicalflow config-urations. These experimental findings were completed by the in-vestigationsofQuirkandKarni [21],providingdetailednumerical resultsofshock–bubbleinteractionproblems.Theyreproducedthe transitionfromregulartoirregularrefraction,shockwavefocusing andtheformationofajettowardsthecenterofthebubble.Fora detailedreviewofSBIseeRanjanetal. [22].

AnewlevelofcomplexitycanbeaddedtothesetupofSBIby replacingtheinertgaswithareactivegasmixture.Astrongshock wave canignitethereactive gasmixturedirectly atthe interface, whereastheadditionalincreaseofpressureandtemperatureinthe shock-focusingpoint isrequiredforignitionatlower shockMach numbers.Twotypeshaveto bedifferentiated:non-premixedand premixedgasmixtures.ReactingSBIofnon-premixedgasmixture was studied by Billet et al. [23]. In their setup a H2 gas bubble

surroundedbyairisshockedtostudytheinfluenceofthevolume viscosity on the bubble evolution and vorticity production. Attal etal. [24] verifiedtheresultsofBilletetal. [23] andfurthermore observed theformation ofa doublediffusionflame inthe bridge region of the shocked bubble. Attal and Ramaprabhu [25] stud-iedsingle-modereactingRMinanon-premixedsetupatdifferent interface thicknesses. They observed shock-induced ignition and mixing enhancement by reshocking the propagating flame. Fur-thermoreshock–flameinteractionincreasesthesurfaceareaofthe flameandtheenergyreleaseandthereforetheburningrate [26]. Massa andJha [27] showedthat smallscales aredamped by the shockwaveandthatthegrowthofRMIandKHIarereduced.

In 2012, Haehn et al. [28] investigated the interaction of a shock wave with a premixed gas bubble, ﬁlled with a stoichio-metricgas mixtureofhydrogen (H2) andoxygen(O2), diluted by

xenon(Xe). Besidestriggering hydrodynamicinstabilities,such as RMI,theshockwave alsoincreasesthetemperatureandpressure, whichinturninducesfasterchemicalreactionratesuptothe igni-tionofthegasmixture.Maximumpressuresandtemperaturesare reachedwhentheshockpassesthebubble.Subsequently,thegas mixturerelaxes andthetwo main parameters controllingthe re-action rate, temperatureandpressure, decrease.The experiments ofHaehn etal. [28] coveredboth ignition typesdeﬂagration and detonation, by varying the shock wave Mach numbers between

Ma=1.34andMa=2.83.

A weak shock wave with Ma=1.34 does not ignite the gas mixture within the experimental timeframe. Compression is not sufficient tostart a self-sustainingchemical reaction. Anincrease oftheshockstrengthresultsinanignitionfollowedbya deflagra-tionwave.ThereactionwavetypechangesforhighershockMach numbers;Haehn etal. [28] observeda detonationwaveforMa= 2.83, even before the shock wave has reached the shock focus-ingpoint.Damköhlernumbersbetween0.25(Ma₌1_.65)and8.00 (Ma=2.83)weredetermined.Haehnetal. [28]concludethatheat conductionplaysanimportantroleatlowershockMachnumbers, and that the Zeldovich mechanism becomes importantat higher shockMachnumbers.Theirconclusionisconsistentwiththetwo limiting cases of shock-induced combustion, the strong and the weak ignition [19].Strong ignition resultsin a detonation essen-tiallyinitiateddirectlybytheshockwave.Weakignitionis charac-terizedbytheappearanceofsmallflamesthatcanundergo transi-tioninto detonationwaves.Severalchemiluminescenceexposures are providedby Haehnetal. [28] todepict thequalitative evolu-tionof thebubbleandreaction processes.Furthermore, quantita-tivedataforthetemporalevolutionofthetransversediameterof thebubble aswellasforthe vortexring diameterarepresented. However, the complex experimental setup implies uncertainties. Haehnetal. [28]estimatetheuncertaintyoftheDamköhler num-ber atthehighestshockMach number(Ma₌2_.83) of upto 50% (Da=8± 4).At the lowestshockMach number(Ma=1.34), 30%

(4)

ofall measurements showednoignition withinthe given experi-mentaltime frame. Hence, numericalstudiesof RSBIcanprovide more certainty and complementary insightinto RSBI phenomena thatcannotbeachievedbypurelyexperimentalwork.

1.4. Scopeofthepresentwork

The presentnumericalstudycomplementsthework ofHaehn et al. [28] and continues the ﬁrst numerical approach to RSBI (Diegelmannetal. [29]).Themainemphasisisplacedonthe gen-eral temporalandspatial evolutionof RSBI, thecomparison with SBI, andthe dependenceofthe bubbleevolution on thereaction wavetype.Inourstudy,theshockMachnumberisvariedbetween

Ma₌2_.13 andMa₌2_.90at a constant initial pressure and tem-perature.Besidesthelimitingcasesofdeﬂagrationanddetonation we studytwo special phenomenaindetail, whichhavenot been discussed before: DDT at Ma₌2_.19 and a double detonation at

Ma=2.50. Haehn et al. [28] observed an effect, which they as-sumeiseitheradoubledetonationorareflectionofmeasurement signals, but the experimental measurement technique did not allowa clearidentification.Ourpresentnumericalstudyconfirms theobservedphysicaleffectandgivesadeeperinsightintothegas composition of thetwo ignition spots duringtheinduction time. Intentionally,wefocuson two-dimensionalconfigurationsasthey facilitate particular analysis and phenomenological investigation. Moreover, in [29] itwas shown that early stages ofRSBI can be wellreproducedbyatwo-dimensionalapproximation.

Thechemicalreactionratesofmostgasmixturesincreasewith pressure.H2–O2reactions,however,showadifferentbehavior [30].

Someintermediatereactionratesareproportionaltothesquareof thepressure,othersarelinearlyproportional [31].Hence,the vari-ationoftheshockMachnumber,ormorepreciselythepost-shock pressure,affectsthechemicalreactionprocessanddeterminesthe occurrenceofeitherdetonationordeﬂagration.

We structurethepaperasfollows: Section2outlines the gov-erningequations,includingmoleculartransportpropertiesfor mul-ticomponent flows and chemical reaction kinetics. Section 3 de-scribesthecomputationaldomainandtheinitialconditionsofour setup. Generalresults arediscussed in Section 4. The spatialand temporal evolution of the RSBI are presented. The effect of dif-ferent types of reaction waves on bubble deformation are com-pared witheach other and withtheir non-reacting counterparts. The chemicalreactionprocess duringshockpassageuntilignition is analyzed in detail. A consistent definition of the dimension-lessDamköhlernumberisusedtoevaluatewhetherhydrodynamic orchemical reactiontime scales dominatetheflow field. Integral quantities,suchasenstrophyorthemolarmixingfraction,are es-timated to assess the effect of the reaction waves on mixingof thebubblegas.In Section5,wediscusstwospecialcasesofRSBI: First the transitionof a deflagration into a detonation wave and second asimulationwithasimultaneousdetonationattwo spots.

Section6presentsacomparisontoexperimentalresultsanda crit-icaldiscussion. Section7summarizesthekeyﬁndings.

2. Numericalmodel

2.1. Navier–Stokesequations

Wesolvethefull setofcompressiblereactingmulticomponent Navier–Stokesequationsinconservativeform

∂ρ

∂

t +

∇

·

(

ρ

u

)

=0 (1)

∂ρ

u

∂

t +

∇

·

(

ρ

uu+p

δ

−

τ

)

=0 (2)

∂

E

∂

t +

∇

· [

(

E+p

)

u]−

∇

·

(

τ

· u− qc − qd

)

− ˙

ω

T=0 (3)

∂ρ

Yi

∂

t +

∇

·

(

ρ

uYi

)

+

∇

· Ji− ˙

ω

i=0 (4)

where

ρ

isthemixturedensityanduthevelocityvector.The iden-titymatrix isgivenby

δ

, totalenergy byE andpressure by p. Yi

arethemassfractionofspeciesi=1,2,...,N,withNbeingthe to-talnumberofspecies.Theheatrelease

ω

˙T andspeciesformation

anddestructionintermsofindividualmassrates

ω

˙irepresentthe

chemicalreactionkinetics.

2.2.Caloricandtransportproperties

Theviscousstresstensor

τ

foraNewtonianﬂuidisgivenby

τ

=2

μ

1 2

∇

u+

(

∇

u

)

T

₋1 3

δ

(

∇

· u

)

, (5)

with

μ

asthemixtureviscosity

μ

= N i=1

μ

iYi/M1i/2 N i=1Yi/M_i1/2 . (6)

Miisdeﬁnedasthemolecularmassofeachspeciesi.The

calcula-tionoftheviscosityofeachspecies

μ

_iisbasedontheChapman– Enskogviscositymodel

μ

i=2.6693· 10−6

MiT

μ,i

σ

i2

, (7)

whereTisthetemperatureand

σ

ithecollisiondiameter.The

col-lisionintegral

_μ,i[32]isdeﬁnedas

μ,i=A

(

Ti∗

)

B+Cexp

(

DTi∗

)

+Eexp

(

FTi∗

)

, (8)

with A₌1_.16145_, B₌₋₀_.14874_, C₌0_.52487_, D₌₋₀_.7732_, E₌

2.16178,F=−2.43787andT_i∗=T/

(

/k

)

i,usingtheLennard–Jones

energyparameter(

/k)iforspeciesi.AccordingtotheFourierlaw,

wedeﬁnetheheatconductionas

qc =−

κ∇

T, (9)

with

κ

asthemixtureheatconductivity,whichiscalculatedfrom

[33]

κ

= N i=1

κ

iYi/M_i1/2 N i=1Yi/M1i/2 , (10)

κ

iis thethermalconductivityofspeciesi.The interspecies

diffu-sionalheatﬂuxq_d[34]isgivenby qd =

N

i=1

hiJi, (11)

withhiastheindividualspeciesenthalpy.ThespeciesdiffusionJi

ismodeledas Ji=−

ρ

Di

∇

Yi− Yi N j=1 Dj

∇

Yj

. (12)

Didescribestheeffectivebinarydiffusioncoeﬃcientofspeciesi

Di=

(

1− Xi

)

_N j=i Xj Di j

−1 , (13)

withXi asthemolefractionofspeciesi. Eq.(13)ensures thatthe

interspeciesdiffusionﬂuxesbalanceto zero.Theconstitutive em-piricallawisused tocompute themassdiffusion coeﬃcientof a binarymixture [33] Di j= 0.0266

D,i j T3/2 p

Mi j

σ

i j2 , (14)

(5)

where Mi j= 2 1 Mi+ 1 Mj and

σ

i j=

σ

i +

σ

j 2 . (15)

Thecollisionintegralfordiffusion

D,ijisgivenby

D,i j =A∗

(

Ti j∗

)

B∗

+C∗exp

(

D∗T_{i j}∗

)

+E∗exp

(

F∗T_{i j}∗

)

+G∗exp

(

H∗T_{i j}∗

)

. (16)

The parameters are deﬁned as A∗=1.06036, B∗=−0.1561, C∗=0.19300, D∗=−0.47635, E∗=1.03587, F∗=−1.52996, G∗= 1_.76474_, H∗₌₋₃_.89411_, and T_{i j}∗₌T_/T_{i j}. T_{i j} have been obtained fromtheLennard–Jonesenergyparametersforspeciesiandjas T_{i j}=

k

i

k

j. (17) 2.3.Equationofstate

Theequationofstateforanidealgasisusedtoclosethe equa-tions

p

(

E,Yi

)

=

(

γ

− 1

)

E. (18)

γ

representstheratioofspeciﬁcheatsofthemixture

γ

= cp cp− R (19) with cp= N i=1 Yicp,i. (20)

The speciﬁc gas constant of the mixture is deﬁned by R₌R_/M_,

withR astheuniversal gasconstant.M isthemolarmassof the mixture M=

_N i=1 Yi Mi

−1 = N i=1 XiMi. (21)

cp,irepresentsthespeciﬁcheatcoeﬃcient

c_p,i=

γ

i

γ

i− 1

Ri. (22)

RiisdeﬁnedasRi=R/Mi.Thetemperatureiscomputedfrom

T= p

R

ρ

. (23)

2.4.Chemicalreactionkinetics

The accurate calculation of chemical reaction kinetics is most important for the precise prediction of combustion effects, such asDDT.The review paper ofOranet al. [19]summarizes several studiesaboutDDT,mainlyoperatingwithone-stepchemical kinet-ics.DDTthroughtheZeldovichgradientmechanismwasobserved, arisingduetothegradientofinductiontimewithin thehotspots in front of the ﬂame, where temperature varies in the range of 600 to 800 K. A precise computationof the induction time and thecorresponding heat releaseistherefore essential foran accu-ratedescriptionofDDT [18]. However, itwasshownthat the in-duction time of detailed mechanisms is larger than for one-step mechanisms [35]andalsolargerthanthetimebetweenﬂame ini-tiationand transitionto detonation,whichrenders numerical re-sultsobtainedwithsimplemechanismsquestionable.Furthermore, importantquantities ofcombustion such as detonation initiation andinduction time in chain-branching kinetics are not correctly reproducedby one-stepmechanisms [36].StudiesofIvanov etal.

[36] revealsigniﬁcantdifferencesbetweenthetemperature gradi-ent that leads to detonation with one-step and detailed mecha-nisms. For the detailed mechanism a much smaller temperature gradient issuﬃcient toignite detonation,whichis inaccordance withthebehaviorofrealcombustiblemixtures [18].

Chemicalreactionkineticsareexpressedbytheheatrelease

ω

˙T

andspeciesformationanddestructionintermsofindividualmass rates

ω

˙i.Thespeciﬁcheatrelease

ω

˙T isdeﬁnedas

˙

ω

T=− N i=1

h0_f,i

ω

˙i, (24) withh0

f,iastheheatofformationofeach speciesi.Massrates

ω

˙i

foreachspeciesareestimatedby

˙

ω

i=Wi NR r=1

ν

ir

r

kf r N j=1 [Xj]ν jr− k br N j=1 [Xj]ν jr

, (25)

withNR asthenumberofreactions,Wi themolecularweight,

r

thethirdbodyeﬃciencyofreactionr,X_jthemolarconcentration, and

ν

_ir and

ν

_irthemolarstoichiometriccoeﬃcientsofthereactant andtheproduct ofreactionr.

ν

ir isthenet stoichiometric

coeﬃ-cient

ν

ir=

ν

ir−

ν

ir. (26)

TheArrhenius lawisusedto calculatetheforwardandbackward reactionrateskfrandkbr.Theforwardreactionratesaredeﬁnedas

kf r=Af rTβf rexp

Ef r RT

, (27)

whereA_fristhepre-exponentialfactor,E_fr istheactivationenergy,

β

fristhetemperatureexponentforeachreactionr [37].The

back-ward reaction rates are calculated by using the equilibrium con-stantsKcr kbr= kf r Kcr. (28) Kcrisgivenby Kcr=

_p◦ RT

νr exp

_S_◦ r,i R −

H_r◦_,i RT

, (29)

withp° asapressureof1atm,

ν

rasthenetchangeinthenumber

ofspeciesinthereaction,

S◦_r,iasthenetchangeinentropy and

H_r◦_,iasthenetchangeinenthalpy.

Furthermore,pressure dependent andduplicated reactions are considered; for this purpose Eq. (27) is modiﬁed. Pressure de-pendence is takeninto account by calculatingtwo forward reac-tion rates kf r₀ and kf r_∞ for the high-pressure and for the

low-pressurelimit,respectively.Ablendingfunctioncomposedofthese high- and low-pressure Arrhenius rate parameters is applied for a smoothpressure dependence.Formoredetails ontheso-called fall-off reactions,the readeris referredto Troe [38].Duplicated re-actionsareconsideredbyextending Eq.(27) to

kf r= 2 i=1 Af riT βf riexp

Ef ri RT

, (30)

Themechanism, whichprovidestheparameters fortheArrhenius law,isessentialfortheaccuracyofthenumericalinvestigationand hastobe chosencarefully.Asshownbytheauthorsinaprevious publication [29], available mechanisms show large discrepancies in ignitiondelay time andpressure sensitivity.A certain number of intermediate reactions, third body eﬃciencies, duplicated and pressure dependent reactions are necessary forthe accurate pre-dictionofthereactionkineticswithinthewiderangeofpressures andtemperaturesconsideredinthisstudy.

(6)

We have chosen the Ó Conaire [39] reaction mechanism for thereactionrateparametersoftheArrheniuslaw.Themechanism is valid for a wide range of pressures (0.05–87 atm) and tem-peratures (298–2700K). 8+N species (two reactants: H2, O2; 5

chain-carrying intermediates: hydrogen radical (H), oxygen radi-cal (O), hydroxyl radical (OH), hydroperoxyl radical (HO2),

hydro-gen peroxid (H2O2);the product: hydrogenoxide (H2O); N inert

gases)and19intermediatereactionsareconsidered,including du-plicated andpressure dependent reactions as well as third-body efficiencies. Third-bodies absorb energy during the two-body re-combinationreactionandstabilizethefinalcombination.The avail-ablemodesforenergystoragecontrol theenergyabsorption.The third-body efficienciesofXearesetidenticaltoargon (Ar),which are provided by Ó Conaire [39]. As the available modes of Ar and Xe are identical, the third-body efficiencies can be assumed to be comparable. Also, the steric factor for monoatomic gases, which accounts for the geometry influence on the collision be-tweenmolecules,issimilar [40].

ThemechanismofÓ Conaire [39]hasbeenusedwidely inthe recent years [41,42]. As part of a preceding validation campaign

[29],theappliedreaction mechanismhasbeencomparedto sim-pler reaction mechanisms. Accurate ignition delay times, crucial for the spatialevolution of the bubble and mixing, can only be achievedwithadetaileddescriptionofchemistrybyasuﬃciently complexreactionmechanism.

2.5. Numericalmethod

The 2nd-order accurate Strang time splitting scheme [43] is usedtosolvethesystemofequations(Eq.(1–4)).TheStrang split-tingschemeseparatesthestiff terms,containingthechemical re-actionkinetics(

ω

˙T and

ω

˙i),fromtheNavier–Stokesequations.This

results in a system of partial differential equations (PDE) andin a stiff system of ordinary differential equations (ODE). We use a finite-volume discretization scheme that applies a flux projec-tion ontolocal characteristicsforthe hyperbolicpart forthePDE system. The Roe matrix required for the projection is calculated forthe full multi-speciessystem [44,45].The numericalfluxes at the cell faces are reconstructed from cell averages by the adap-tivecentral-upwind6th-orderweightedessentiallynon-oscillatory (WENO-CU6)scheme [46].Theschemeusesanon-dissipative6th -ordercentralstencilinsmoothflowregionsandanon-linear con-vex combination of 3rd-orderstencils in regions with steep gra-dients. Time integration is performed by the 3rd-order strongly stable Runge–Kutta scheme,developed by Gottlieb and Shu [47]. Ournumericalmodelhasbeentestedandvalidatedforshock in-ducedturbulentmulti-speciesmixingproblemsatfiniteReynolds numbers [13,48–50] and for shock–bubble interactions including chemistry [29].The stiff ODE,governing the specific heatrelease andmassratesforeachspecies,isseparatelysolvedbya variable-coefficient ODE solver using 5th-order backward differentiation formulas [51].

3. Computationalsetup

We study RSBI within a two-dimensional rectangular domain withasymmetryplaneatthecenteraxisofthebubble,see Fig.1. Inﬂowboundaryconditionsareimposedattheleftdomain bound-aryand outﬂowboundary conditionsat the rightandupper do-mainboundaries.Thedomainsizeissetto32.5r_{× 10.5}r,withr

astheinitial bubbleradius.Thedistancebetweenthebubbleand domainboundariesare chosen sufficientlylarge to avoidartifacts duetoshockreflections.TheCartesian gridintheregionof inter-estisrefinedby afactorof25compared tothecoarseouter grid toreducecomputationalcosts.Adetailedgridstudycanbefound inourpreviouspaper.Wedemonstratedgridconvergenceby com-paringfourdifferentgridresolutions,withcellsizesof

xy=234,

117,59and29μminthehighresolutionpart.Thesimulationsare performedataCFL-numberof0.3.

Thegas bubbleisﬁlled withH2,O2 andXe ina

stoichiomet-riccompositionof2/1/3.67massfractionsandsurroundedbypure Nitrogen(N2).ThebubblediameterissettoD=2r=0.04m.The

heavy inertgasXe isused to increase thedensityof thebubble, leadingto an Atwood numberof A=0.476. The gascomposition of our domain and the bubble diameterare identical to the ex-perimentalsetup ofHaehnet al. [28].A sharpandfullyresolved interfacebetweenthebubblegasanditssurroundingisdeﬁnedin termsofthemolarfractionofN2

XN2=

tanh

((

x2₊_y2_{− r}

)

ξ

)

₊₁

2 , (31)

withr as the radius of the bubble and

ξ

asparameter for con-trolling steepness, which is set to

ξ

₌ 20,000. The molar frac-tion(X=1− XN₂)insidethebubbleisdistributedamongthethree

gases,ensuringthestoichiometricmixturewitharelative compo-sitionof2/1/3.67(H2/O2/Xe).

The shock wave is initialized on the left side of the bubble. The pre-shockstate is deﬁnedby T0=350K and p0=0.50atm.

The shockMach number is varied between Ma=2.13 and Ma= 2.90. The post-shock thermodynamicsstate is given by standard Rankine–Hugoniotconditions

ρ

N2=

ρ

N2

(

γ

N2+1

)

Ma 2 2+

(

γ

N2− 1

)

Ma2 , (32) uN2=MacN2

1−

ρN

2

ρ

N2

, (33) p_N₂=p0

1+2

γ

N2

γ

N2+1

(

Ma2_{− 1}

₎

, (34) where cN2=

_γ

N2p0/

ρ

N2. Variables indicating post-shock

condi-tionsaremarkedwithaprime.

Note that the initial parameters of our setup slightly deviate from the experimental pressure and temperature. To avoid that

(7)

thepressurepeakofthedetonationfrontisoutsideofthe valida-tionrangeofcurrentlyavailablereactionmechanismsfordetailed H2–O2 reaction kinetics,we slightlydecrease the initial pressure

andincrease the initial temperature as compared to the experi-ment to achieve a similar reaction behavior. We believe that it isimportantforfurthernumericalinvestigations tooperateinside thevalidatedrangeofthereactionmechanism.

4. Resultsanddiscussion

The present numerical investigation of RSBI covers different reaction wave types triggered by shock Mach numbers between

Ma₌2_.13 andMa₌2_.90. Deﬂagration is induced by the lowest shockMachnumberofMa1=2.13.Increasingshockstrengthleads

tothree differenttypesofsupersonicreaction waves:Ma2=2.19

induces a deﬂagration, which transitions to a detonation. Ma3=

2.30 and Ma5=2.90 immediately cause detonations behind the

shockwaveatthedownstreamorupstreampoleofthebubble, re-spectively.AshockMachnumberofMa4=2.50leadsto anearly

simultaneous double detonation in two bubble regions. Tempo-ral and spatial bubble evolution, enstrophy production and mix-ing are strongly affected by the reaction wave type, which we discusscomprehensivelyinthefollowingsections.Thesimulation witha shock Mach number of Ma2=2.19 is excluded from the

discussion,as the globalbubble dynamicsare nearly identicalto

RSBI at Ma3=2.30. The transition process will be discussed in Section 5.1andthe doubledetonation atMa4=2.50 willbe

dis-cussedin Section5.2.

4.1. Globalbubbledynamics

The qualitative inﬂuence of the chemical reaction kinetics on thetemporalevolutionofSBIisshownin Fig.2.Thecontourplots of thetemperature inside the bubbleshow the compressionand propagationofthereactionfront.Theupperpartofeachsequence showsthereactingsimulation,thelowerpartprovidesresultsfor thenon-reactingsimulationatthesameshockMachnumber.For clarity we ﬁrst compare simulations with shock Mach numbers

Ma1=2.13,Ma3=2.30andMa5=2.90.

We refer to the lower part of the ﬁrst sequence of contour plotsforageneraldescriptionofthecharacteristicstagesofbubble evolution.At t=50μs andt=86μs,the shockwave propagates throughthebubbleandcompressesthegasmixtureinside.Before theshockwavehaspassed,ﬁrstinstabilitiesontheinterface start to arise (t=86μs). At t=120μs, the roll-up of the bubble has started,primaryvorticesformandsecondaryinstabilitiesgrowdue toshearatthematerialinterface.Finallyatt=500μs,thebubble gasshowsa highdegree ofmixingwiththe surroundinggasN2.

Thetwomainvorticesarefullydevelopedandconnectedoverthe bridgeregionattheupstreampoleofthebubble.Secondaryvortex

Fig. 2. Temperature contour plots of SBI: upper parts show reacting SBI, lower parts the inert counterparts. Ma 1 = 2 . 13 induces a deﬂagration wave, M a 3 = 2 . 30 and M a 5 = 2 . 90 a detonation wave.

(8)

structures are clearly visible at the outer interface. The general evolutionofinertSBIandthecharacteristicstagesaresimilar.The differentshockstrengthsleadstoarangeofpropagationvelocities of the shock waves,which in turn shift the overall bubble evo-lution intime.Furthermore, ahighershockMach numbercauses ﬁnerstructuresinthelong-termdevelopmentoftheSBI.

Similarity of evolution at different shock Mach numbers can-not be observed when chemical reaction kinetics are taken into account.Asthereactionsarestronglypressuresensitive,theshock Machnumberaffectstheinductiontimeandthesubsequent reac-tionprocess.Theﬁrstandweakestshockwave (Ma1=2.13)leads

todeﬂagration.Thegasmixtureisignitedshortlybeforet=120μs, as shown in the upper row of Fig. 2. The reaction front propa-gates slowlythrough the bubble gas; even at t=500 μs the re-action front hasnot yet reachedthe upstream pole.Thus bubble structures ofreacting andnon-reactingSBI arestill similar, espe-ciallytheouterinterfaceswithevolvingKelvin–Helmholtz instabil-itiesshow thesamecharacteristics.Merely,theoverallbubble ex-pansionincreasesduetoisobaricheatingoverthereactionfront.

Whenwe increasetheshockMachnumbertoMa3=2.30,the

bubble dynamic changes distinctly. The reactive gas mixture ig-nites earlier,followedby a detonationwave, depictedinthe sec-ondrowof Fig.2.Thesupersonicreactionwavepropagateswithin approximately

t=10μsthroughthebubble.Strongheat release and densitydecreaseresult ina rapid andsigniﬁcant bubble ex-pansion. When the detonation wave reaches the interface of the bubble, vorticity is generated, with the opposite sign compared to the vorticity induced by the initial shock wave. As a conse-quence, the growth of the secondary instabilities is decelerated, which has a signiﬁcant effect on mixing. The contour plots at

t=120μsandt=500μsshowthebubbleafterthereactionwave haspropagatedthroughthereactivegasmixture.Comparisonwith the inert counterpart reveals different growth rates and charac-teristics.Furthermore,thedetonationwaveampliﬁestheN2-jetat

thesymmetryplaneatthedownstreampoleofthebubble.Hence the bridge region, connecting the two main vortices, vanishes completely.

A further increase of the shock Mach number to Ma5=2.90

shortens the inductionsignificantly. The strong shockignites the gasmixturedirectlyattheupstream poleofthe bubble,followed by a detonationwave. Thedetonation wave mergeswiththe ini-tial shock wave inside the bubble and subsequently propagates through the unshocked gas mixture. Comparison with the inert counterpartshowsthatthedetonationwavepropagatesmorethan twice as fast asthe initial shock wave, which significantly influ-encesthespatialbubbleevolution.Similarlytothesimulationwith

Ma3=2.30, some of the secondary instabilities are suppressed;

however,theN2-jetislessampliﬁedandthebridgeregionis

pre-served evenin thelong-termevolution.The bubbleispenetrated byasingledetonationwave,whereasthelowershockMach num-ber of Ma3=2.30 induces a detonation wave that propagatesin

upstream directionthroughthe pre-shockgasmixture.The inter-face differs from the other simulations, the KHI are not aligned along theouterinterface.Furthermore,thetimestepatt=120μs showssecondaryRMI arisingfromtheKHI. Differentevolutionof primaryandsecondaryinstabilities,thebridgeregion,andthe spa-tialbubble expansionhavea signiﬁcant effectonintegral quanti-ties,suchasmixingandenstrophyproduction.

The vortex Reynolds number Re₌

0/

ν

, deﬁned by Glezer [52],where

0 istheinitialdepositionofvorticityand

ν

the

kine-maticviscosityoftheinterface,amountsto1.4× 105_for_Ma₌₂_.₁₃

and increases with higher shock Mach numbers up to 2_.4_{× 10}5

(Ma=2.90).ThecriticalReynoldsnumberfortransitionfrom lam-inartoturbulentﬂowis104_to₂_{× 10}4_(Dimotakis_[53]_)._All

conﬁg-urationsexceedthismixing-transitionReynoldsnumberbyatleast oneorderofmagnitude.

Fig. 3. Normalized transverse bubble diameter for different shock Mach numbers. — : reaction; – : no reaction; : Ma 1 = 2 . 13 , : Ma 3 = 2 . 30 , : Ma 5 = 2 . 90 .

4.2.Transversebubblediameter

Thetransversebubblediameter

y=

y/D0normalizedbythe

initialbubblediameterD0isusedtomeasuretheimpactof

chem-icalreactionprocesseson thelarge-scaleevolution ofthebubble. Thebubble diameter

y ismeasured based ona thresholdvalue

ofthexenonmassfractionofYXe=0.01.

Figure3showsthe temporalevolutionof

y fortheinertand

the reacting simulations at three different shockMach numbers. Fortheinertsimulations(dashed-dottedlines)weobserveanearly linear increase in the bubble diameter. Some variation can be found inthelong-term evolution:the weaker theshockstrength the smaller the streamwise expansion, which leads to a slightly largertransversebubblediameter.AtthehighestshockMach num-berofMa5=2.90,theevolution levels outmuch earlierthanfor

theother shockMach numbers.Note that theroll-up of the pri-maryvorticesleadstoa wave-liketemporalgrowthof

y.Inthe

long-termevolutionthemainvorticesarefullydevelopedandthe bubblegasrotatesaroundthevortexcores,whichresultsina ﬂat-teningofthetransversebubblediameterevolution.

The deﬂagration wave triggered at a shock Mach number of

Ma1=2.13affectsthenormalizedtransversebubblediameteronly

with respect to the long-term evolution. The propagation veloc-ityofthe reactionfront is low comparedto the evolutionof the hydrodynamicinstabilities.Thedensityincrease overthe reaction front accompaniedby a spatialexpansion leadstoa slight diver-gencefromtheinertcounterpartaftert₌300μs.Simulationswith shock Mach numbers of Ma3=2.30 and Ma5=2.90 exhibit

su-personicdetonationwaves,whichhaveasigniﬁcantly stronger ef-fect onthe normalizedtransverse bubblediameter. Therapid re-actionleadstoaninstantaneousexpansionofthereactedgasand a sudden increase of

y up to 175% of the initial bubble

diam-eter.The detonationwave atMa3=2.30propagates inupstream

andorthogonal direction ofthe ﬂow ﬁeld. Thus, an overshoot in thetransversebubblediameterisvisibleatt=180μs,which de-creasesovertimetoaslightlylowervalue.ThehighershockMach numberofMa5=2.90showsanearlierincreaseofthetransverse

bubble diameter and levels out at a higher value compared to the simulation with a shock Mach number of Ma3=2.30. The

higher compression leads to a higher temperature and therefore toalargerspatialexpansionofthebubblegas.

4.3.Identiﬁcationofthereactionwavetype

Tooutlinethedifferentfeaturesofdeﬂagrationanddetonation wavesweanalyzetheevolutionofcharacteristicparametersacross

(9)

Fig. 4. Pressure and gas composition across the fully developed reaction wave front. — : pressure; – : Y H 2O ; — : Y H O 2 .

thereactionfront. Figure4showsthepressureandthemass frac-tionsoftheradicalHO2 andtheproductgasH2Oacrossthe

reac-tionwavesforMa1=2.13 (deﬂagration)andMa3=2.30

(detona-tion)atatimeandposition,whenthereactionwaveisfully devel-opedandtheinitialshockwavehaspassedthebubble.(t=90μs forthedetonationwaveandt₌300μsforthedeﬂagrationwave.)

Figure4(a)showsdataforthedeﬂagrationreactionwave.The reactionfront propagatesfromrightto left anditslocation coin-cides withthe peaks of the HO2 mass fraction. Accompaniedby

thepeakarapidincreaseoftheproductgaseousH2Oisapparent.

Thepressureisnotaffectedbythechemicalreactionandremains constantacrossthereactionfront.

A different evolution is observed for the supersonic reaction wave at Ma3=2.30, see Fig. 4(b).In addition to the increase of

theproductgaseousH2OandtheintermediatespeciesHO2across

thereaction front, alsothe pressure exhibitsa pronounced peak, whichiscausedbytheshockwaveprecedingthedetonationwave. Thepressuredecreasesbehindtheshockwavebutlevelsout ata largervaluecomparedtothedeﬂagrationwave.Moreover,the am-plitudesof theHO2 mass fractionpeaksinthereactionzone

dif-fer. The detonation wave shows a signiﬁcantly higher amount of HO2andabreakdownacrossthereactionzoneindicatingthatthe

thirdexplosionlimit iscrossed andhigh-pressurereactions dom-inate [54].Bothreactionwaves resultin thesameamountof the productgasH2O.

Figure5showsthemassfractionofH2O2 forthetwodifferent

reactionwavetypes.Thespatialcoordinate

ξ

denotesthedistance fromthe reactionfront, which propagatesfromright toleft. The peakofYH₂O₂ indicates areactionabovethethirdexplosionlimit.

Fig. 5. H 2 O 2 mass fraction across the fully developed reaction front. — : Ma 1 = 2 . 13 (deﬂagration) and — Ma 3 = 2 . 30 (detonation).

HO2 collides with H2 forming eitherH2O2 ordirectly H2O. [54].

The strongreduction ofYH2O2 on theright sideof theplot

iden-tiﬁes the bubbleinterface and consequently the boundaryof the reactionzone.

4.4. Damköhlernumber

Theflow field ofRSBIis affectedbyhydrodynamiceffectsand chemicalreactionkinetics.TheDamköhlernumber,definedasthe ratioofthehydrodynamicandchemicalreactiontimescales, Da=

τ

h

τ

r,

(35)

indicateswhicheffectdominates.Da>1characterizesaflowfield mainlydrivenbychemicalreactions,Da<1impliesadomination of the hydrodynamiceffects. The two time scales are definedas follows:

τ

h= 1

|

ω

|

, (36)

τ

r=

τ

ign+ D0 2VRW. (37)

Thecharacteristichydrodynamictimescale

τ

hisdeﬁnedbythe

to-talvorticity

ω

,averagedfromtheﬁrst contactoftheshockwave withthe bubble until the reaction wave has propagated through the bubble. The chemical reaction time scale

τ

r consists of two

timeintervals:

τ

ignistheperiodfromtheﬁrstcontactoftheshock

withthe bubbleuntil ignition,andD0/(2VRW) isthe time the

re-action wave needs to propagate through half of the initial bub-bleshapewithD0 astheinitialbubblediameter.Thepropagation

velocity of the deﬂagration wave is sensitive to temperatureand pressure,considered by apowerlaw expressionasintroduced by Dehoe [55] VRW =SL0

_T T0

β1

_p p0

β2 , (38)

withSL0 denotingthelaminarburningvelocityatreference

condi-tions (T0 and p0), available inrecent literature [55]. Temperature

Tandpressurepare takenfromthe hotspotshortlybefore igni-tion.Theparameters

β

1 and

β

2 are1.54and0.43, respectively,as

we aredealingwitha stoichiometricmixture [56]. Wealso com-puted the propagation velocity directly from the simulation and found good agreement tothe literature. As the varying tempera-tureandpressuredistributionsinsidethebubbleleadtoarangeof differentpropagationvelocities,wedecidedtousethevelocity cal-culatedfrom Eq.(38).Theobtainedvelocityprovedtobea reason-able estimatedvalue forthe propagationvelocity. Fordetonation the propagationvelocity ofthe reactionwave is morestableand thereforedetermineddirectlyfromthesimulation. Table1provides

(10)

Table 1

Damköhler numbers and characteristic time scales for different shock Mach numbers. Ma [ −] τh [ s ] τr [ s ] Da [ −] 2.13 4 . 926 × 10 −4 1 . 739 × 10 −3 0.283 2.19 1 . 227 × 10 −4 9 . 291 × 10 −5 1.321 2.30 1 . 094 × 10 −4 8 . 278 × 10 −5 1.322 2.50 1 . 286 × 10 −4 7 . 499 × 10 −5 1.715 2.90 3 . 086 × 10 −4 4 . 499 × 10 −5 6.859

theDamköhlernumbersforthedifferentinitialshockMach num-bers,includingtheonesforMa2=2.19andMa4=2.50,whichare

discussedin Section5.

ThedeﬂagrationwaveinducedatalowshockMachnumberof

Ma1=2.13leads toa Damköhler numberofDa=0.283. The

hy-drodynamictime scaledominatestheflow field,bubbleevolution andthegrowthofsecondaryinstabilitiesaremainlydrivenby hy-drodynamic effects.Asshownin Fig.2,theinfluenceofthe reac-tionfrontisminor.

AthighershockMachnumbers,chemicalreactionsstarttoplay a crucial role for the overall bubbledynamics. The fast propaga-tionvelocityofthedetonationwaveshortensthechemicalreaction timescale

τ

randleadstoaDamköhlernumberofDa=1.322for

a shock Mach number of Ma3=2.30. The primary vortex region

andthe growth ofsecondary instabilities are directlyaffected by thedetonationwave.AtthehighestshockMachnumberofMa5=

2.90,evolution isentirelydominatedby thedetonationwave. Af-ter it merges withthe initial shockwave, it determines the spa-tialevolutionofthegasbubble.RMI,theprimaryvorticesandthe secondaryinstabilitiesemerge undertheinfluenceofthereaction wave.Theignitionattheupstreampoleofthebubbleshortensthe reactiontimescale.Thesinglereactionwavereducesthevorticity productioncomparedto thedetonationwave,which originatesat thedownstreampoleofthebubbleandthereforeincreasesthe hy-drodynamictimescale.Thereactionwavedominatestheflowfield, which finds expression in a significant increase of the Damköh-ler numberup to Da=6.859.At a lower shockMach number of

Ma3=2.30,the earlyRMI evolvesin theunburnt gasmixture of

thebubble,which leadsto alower Damköhlernumbercompared tothesimulationatMa5=2.90.

4.5. Enstrophygeneration

Weusetheenstrophy

ε

=

S

ω

2_d_x_d_y ₍₃₉₎

to determinetheinfluenceofthe differentreactionwavesonthe vorticityproduction. Figure6outlines theenstrophy overtime for reactingandinertSBI.Theenstrophyiszerountiltheshockwave reachestheupstream poleofthebubble.Baroclinic vorticity pro-ductionleadstoanincreaseduringtheshockwavepassage.Afirst localmaximuminenstrophyisreachedaftertheshockhaspassed half ofthe bubble,an effectthat can be observedfor all simula-tions. Thereafter,aslightdecayisvisible,followedbyanother in-creaseduetoshockfocusingandshockreflectionsattheinterface. Theenstrophygraduallydecaysafterthepassageoftheshock.The same patternis observed forall inertsimulations, independently oftheshockMach number;onlyoverallenstrophy levelsdifferas strongershockwavesgeneratemoreenstrophy.

The deﬂagration wave induced by a shock Mach number of

Ma1=2.13hasnonoticeableinﬂuenceontheenstrophy,the

vari-ationbetweenthereactingandinertsimulationsisnegligible.The detonationwavesinducedbyashockMachnumberofMa3=2.30

producesigniﬁcantamountsofadditionalvorticity,whichleadsto

Fig. 6. Enstrophy. — : reaction; — : no reaction; : Ma 1 = 2 . 13 , : Ma 3 = 2 . 30 , : Ma 5 = 2 . 90 .

adistinct enstrophy peak, see Fig.6.The elevated enstrophy lev-els persist for about 50 μs. The highest shock Mach number of

Ma5=2.90showsa differentbehavior.As theignition spotis

lo-catedatthe upstream poleofthe bubbleand asthemixture ig-nitesimmediatelyaftertheﬁrstcontactoftheshockwave, enstro-phyproduction is dominated by the detonationwave. Therefore, we have two enstrophy peaks of similar magnitude, one during the shockwave passageof the upstream part of thebubble and one duringthe passageof the downstream part.Thereafter, sim-ilar to thesimulation at Ma3=2.30, theenstrophy decays faster

comparedtotheir inertcounterparts.The inertSBI,showninthe contour plots in Fig. 2, is characterized by several smallvortices attheouter interface,eveninthelong-termevolution.The react-ingSBIshowsamuchsmootherﬂowﬁeldwithfewervortices.The detonation waves of reacting SBIsdecelerate the growth of sec-ondaryinstabilitiesandreducestheappearanceofsmallervortices asit inducesvorticity withoppositesigncompared tothe vortic-ityproducedbytheinitialshockwave.Furthermoretheincreased diffusionacross thereaction frontdamps thegrowthrateof sec-ondaryinstabilities.The shockwave MachnumberofMa5=2.90

revealsanadditionaleffect:enstrophy productionduringthe first halfoftheshockwavepassageforthereactingSBIishigherthan fortheinertSBI.Partsofthedetonationwave arereflected,when it merges with the initial shock wave. The reflected wave pro-ducesadditionalvorticity atthe internalinterfaceinside the bub-ble.Duringthesecondhalfoftheshockwavepassage,theinertSBI showsa largerenstrophy production.The densitygradient across theshockwave ishigherthanthe gradientacross thedetonation wave, leading to a higher enstrophy production. The vorticity of thereflected shockwave ina reactingSBI, increasing the enstro-phyduringthe first partof the shockwave passage, hasalready decayedatthisstageofSBI.

4.6.MixinginRSBI

The shock–bubble interaction provides a complex flow field, where RMI and KHI induce local spots of high mixing rates. Tomkins et al. [57] identified three main regions of mixing: the mainvortices,theouterinterfaceincludingKHIandthebridge re-gion,whichconnectsthetwomainvortices.Thelattercontributes up to 40% to the mixing. Toestimate theimpact of the reaction wavesonthemixing,weusethemolecularmixingfraction(MMF), definedbyDanckwerts [58]as

(

t

)

= _∞ −∞

XN2XXe

dx _∞ −∞

XN2

XXe

dx . (40)

(11)

Fig. 7. Molecular mixing fraction. — : reaction; — : no reaction; : Ma 1 = 2 . 13 ,  : Ma 3 = 2 . 30 , : Ma 5 = 2 . 90 .

Themolarmixingfractioncanbeinterpretedastheratioof molec-ularmixingto large-scale entrainmentby convective motion.We plotthetemporal evolutionof

(t) forreactingandinert simula-tionsin Fig.7.

The inertsimulations showa lineargrowthofthemolar mix-ingfraction,theslopeincreaseswithhighershockMach number. Anincreaseinshockstrengthleadstohigherenstrophyproduction andfastergrowth ofsecondary instabilities, whichenhance mix-ing.The reacting counterparts show a different behavior. In gen-eral, the mixing is reduced by the reaction waves, independent oftheir type.The deﬂagration wave inducedby thelowest shock Machnumberleadstoadecreaseofup to30%.Mixingisaffected

after approximately t=270 μs, when the reaction wave reaches partsof theinterface andthemain vortices. However,the bridge region remains unaffected,asthe propagationvelocity ofthe de-ﬂagrationwaveistoolow.

Thedetonationwaveaffectsallthreemainmixingregions.The MMFisreducedbyupto50%forMa3=2.30aswellasforMa5=

2_.90. Besides the reduction of mixing in the region of the main vorticesandattheinterface,whicharealreadyaffectedbya defla-gration waveatalower shockMach number,thebridgeregionis alsoinfluencedbythedetonationwave.Thedetonationwaves de-celerate thegrowthof secondaryinstabilities. Especiallythe bub-bleevolutionatthehighestshockMachnumber,see Fig.2,shows dampingoffinestructures,whichexplainsthehigherreductionof theMMF.

Figure 8outlines the mixingprogress inthe long-term evolu-tion at t=500 μs for three different shock Mach numbers. The inertSBIsshowalreadyahighdegree ofmixing,whereas the re-actingSBIsarecharacterizedbylargeareasofunmixedbubblegas. TheinertSBIsshowhighermixingforincreasingshockMach num-ber, whichisinaccordanceto theMMF plottedin Fig.7.The re-actingSBIsfollowasimilartrend:ForlowshockMachnumbersof

Ma=2.13(Fig.8(a)) andMa=2.30 (Fig. 8(b)) we observelarger areas of unmixed bubble gas, whereas the highest shock Mach numberofMa=2.90(Fig.8(c))showsahigherdegreeofmixing.

5. Specialcases

Two simulations of RSBI contain hydrodynamic and chemical features that have to be discussed in detail. We observed DDT for Ma2=2.19 and a double detonation with different reaction

branchesforMa4=2.50.

Fig. 8. Mixing of N 2 and Xe for inert and reacting SBI at t = 500 μs for different shock Mach numbers. Upper parts show the reacting SBI, lower parts the inert counterpart.

(12)

Fig. 10. Detailed contour plot of RSBI at Ma 2 = 2 . 19 shortly after DDT: upper part shows the temperature, lower part the pressure. The arrow shows the line slice for the data of Fig. 11 .

5.1. Deﬂagration-to-detonationtransition

(

Ma2=2.19

)

Figure 9 shows contour plots forMa2=2.19 duringthe early

stageofignition.Temperature(upperpartoftheplot)andpressure (lowerpart)areparameters thatillustratetheignitionand transi-tionprocess.Att=93.0μstheshockwavestillpropagatesthrough thebubblegas,andinstabilitiesattheinterfacestarttogrow.The gasmixtureignitesnearthedownstreampoleofthebubbleat ap-proximately t=93.2 μs.A subsonic deﬂagration wave propagates through thebubblegasuntilt₌94_.9μs.At t₌94_.9μs the tran-sitionintoadetonationwavestartsinthelower regionofthe re-action front.The remaining contourplots showsa fastgrowthof the supersonicdetonationwave out ofthe deﬂagration front and acharacteristicsteeppressurepeakacrossthedetonationreaction wave.

Detailed contour plots to illustrate how the detonation wave evolvesfromthedeﬂagrationwaveareoutlinedin Fig.10.The up-per partof theﬁgure showsthe temperatureandthe lowerpart thepressure.Atemperatureisoline(T=2000K)visualizesthe re-actionfront.

The temporal evolution of the characteristic thermodynamic propertiesduringthetransitionisplottedin Fig.11.Theplotsshow the variationsof pressure,theradical Hconcentration (Fig.11(a)) and the temperature (Fig. 11(b)) for seven timesteps from t= 94.6 μs until t=96.0 μs. The coordinate

ξ

is obtained by a ro-tationof theoriginal coordinatesystem. As a resultofthe

trans-Fig. 12. Reaction front marked by black temperature isocontour, shortly be- fore transition to detonation. High-pressure region denoted by red isocontour ( p > 5 bar).

formation

ξ

coincides with the propagation direction of the re-actionwave. The reaction wave propagatesfrom the rightto the left. We observed an increase ofthe pressure peak accompanied byadecreaseofHinthereactionzone,whichischaracteristicfor DDT.The Hradical appears within thedeflagration wave, follow-ingthe detonationwave.The peakofHcan befound behindthe shock wave, after the sudden increase in temperatureand pres-sure.Furthermore,DDTcanbeidentifiedbythesteepeningofthe temperatureprofile,its peak atthereaction frontandthe higher temperatureoftheproductgas.Thesefindingsareconsistentwith observationsofIvanovetal. [36]andLibermanetal. [15,18].

The deﬂagration front propagates in semicircular direction throughthebubblegas,hencethequestionariseswhythe transi-tiontodetonationoccursatspeciﬁcareasofthereactionfront.The reasoncanbefoundbythedetailedanalysisofthepressure distri-butioninfrontofthereactionfront,whereDDToccurs. Figure 12

showsthe reaction front shortly before thetransition to detona-tion.Thereactionzone isindicated bytwo blackisocontours,the red isocontourdepictsthe region of highpressure. Dueto shock focusingof theinitial shockwave atthe downstreampole ofthe bubble,aregion ofhighpressure(p> 5bar) exists.The reaction front propagatesinto thisregion and high-pressurereactions are promoted,whichsupportthetransitiontodetonation.

Asthedeﬂagration wavepersists onlyfora fewmicroseconds, theoverallbubbleevolutionoftheRSBIatMa2=2.19issimilarto

thatatashockwaveMachnumberofMa3=2.30.Thenormalized

transversebubblediameter,theenstrophyproductionandthe mo-larmixing fraction are nearly identical. The Damköhler numbers

Fig. 11. Pressure, temperature and radical H concentration in the reaction front during the transition to detonation for seven conservative time steps with t = 0 . 2 μs starting at t = 94 . 6 μs.

(13)

Fig. 13. Contour plots of RSBI at Ma 4 = 2 . 50 : upper parts shows the reacting SBI, lower parts the inert SBI.

amount to Da=1.321 (Ma2=2.19) and Da=1.322 (Ma3=2.30)

indicating a chemically dominated ﬂow ﬁeld. To ensure that the observedDDTisnot causedby numericalartifacts,theRSBIwith ashockMachnumberofMa2=2.19wasrepeatedatacoarsened

andreﬁnedgridresolution.Weusedthesamegridsthat were al-readyappliedforthegridconvergencestudyinourpreviouspaper

[29]withcellsizesof

xy = 234,117 and59μm.Allsimulations

showthesameinductiontime,ignitionspotandlocationofDDT.

5.2.Doubledetonation

(

Ma4=2.50

)

The shockwave of Ma4=2.50 induces two detonationwaves.

Oneoriginatesnearthedownstreampoleofthebubbleandoneat theupstreampole. Figure13showstheignitionspots,the propa-gationandinteractionof thedetonation waves.The tworeaction wavespropagatetowards eachother,which leadstoa rapid con-sumptionofthereactivebubblegas.Thereactiontimescaleis sig-nificantlyshortened,leadingtoanincreaseoftheDamköhler num-bertoDa=1.715,see Table 1.Specific conditionshavetobe sat-isfiedtocauseadoubledetonation:Theignitiondelaytimeinthe first reaction spotat the upstream pole of the bubblehas to be longerthan in thesecond spotnearthe downstream pole. How-ever, the ignition has to be nearly simultaneous in the absolute timeframe.

The reactionregionattheupstreampoleischaracterized bya slowincrease of theintermediate gas mass fractions directly be-hind the shock wave, beginning at t=28 μs. After an induction timeof about47μs the gasmixtureignites andformsa detona-tionwave.The secondignitionspotshowsadifferentbehavior;a strongreactionisdirectlyinducedaftertheshockwavehaspassed att=68μs.Theinductiontimeonlyamountstot=6μs,whichis muchshorterthanfortheﬁrstignitionspot.Thehigher compres-sionandtemperatureatthedownstream poleofthebubblelead toa fasterignition. In the absolutetimeframe both spots trigger ignition,followedbyadetonationwaveatapproximatelythesame time,whichleadstoadoubledetonation.

Figure 14outlines the temporal evolution ofthe intermediate speciesH, O, OH during the induction time in the two reaction zones. The solid lines denote the ignition spot at the upstream poleofthebubble,thedashedlinesshowdataforthedownstream pole.As the chemical reactions are highlypressure sensitive the highertemperatureandpressureatthedownstreampoleleadsto afasterproductionoftheradicals andtoa shorterignitiondelay time,comparedtotheupstreampole.

6. Discussion

6.1.Transversebubblediameter

OursetupisbasedontheexperimentalinvestigationofHaehn etal. [28].Theyobserved,similarlytoourwork,differentreaction wavetypesbyvaryingtheshockMachnumberbetweenMa=1.34

Fig. 14. Maximum mass fraction of intermediate products inside the two ignition spots. = H , = O , = OH . — : 1st ignition spot at the upstream pole; — : 2nd ignition spot at the downstream pole.

and Ma=2.83. In the following section, we compare the results of the two studies. We are well aware that a comparison of two-dimensionalsimulations withtheexperimentofHaehn etal.

[28]isonlyreasonablealongearlystagesofevolution [29].Haehn etal. [28] providetheDamköhlernumber,thetransverse normal-izedbubbleandthemainvortexdiameterforseveralexperimental setups.Wewillusethetransversebubblediametertocomparethe generalevolution ofthebubbleandthe inﬂuenceof thereaction waves.Theexperimentsexhibitdeﬂagration forshockMach num-bersofMa₌1_.63andMa₌2_.07anddetonationforashockMach numberofMa=2.83.Haehnetal. [40]alsoobserveddetonations forshockMachnumbers smallerthan Ma=2.83, however, with-out providing quantitative data. We compare the experimental results for Maexp=2.07 and Maexp=2.83 with our results

ob-tained at Manum=2.13 and Manum=2.90. Figure 15 shows the

evolutionofthenormalizedtransversebubblediameterforaRSBI witheitheradeﬂagration wave (a)oradetonationwave (b).The normalizedtimet∗followsthedeﬁnitionofHaehnetal. [40] t∗=

_τ

t

n,

(41)

wheret isthetime measured fromtheﬁrstcontactoftheshock wave with the bubble.

τ

n is deﬁned as D0/Wi, with Wi as the

incidentshockwavespeed.

Asdiscussedin Section4.2,theslowlypropagatingdeﬂagration wave leadsonly toaslightincrease ofthebubblediameter.Inert andreactingSBIshowasimilarevolution,conﬁrmedbythe exper-imentalaswellasthenumericalresultsin Fig.15(a).The numer-ical data show a larger normalizedbubble diameter inthe long-termevolution,whichisattributedtotwo-dimensionaleffects.

Also the data for SBIs, which induce a detonation wave, are invery goodagreement. Thepropagationvelocity ofthe reaction wave, the spatial expansion of the bubble and the peak of the

(14)

Fig. 15. Normalized transverse bubble diameter. Comparison between simulations and experimental data of Haehn et al. [28] . Reaction: — and ; No reaction: — and . (a) Deﬂagration Ma = 2 . 13 / 2 . 07 [28] , (b) Detonation Ma = 2 . 90 / 2 . 83 [28] .

transverse bubble diameter are nearly identical. Figure 15(b) re-vealsagoodmatchofthesteepbubblediameterincrease, indicat-ingthatthepropagationofthereactionwaveinsidethegasbubble is well reproduced by thetwo-dimensional simulations. Similarly to the deﬂagration setup, the bubble expansion differs from the experimentalresultsinthelong-termevolution.

6.2. Doubledetonation

As mentionedbefore,Haehnetal. [28] alsoperformed experi-mentsatintermediateshockstrengths.Foronespecificsetup,their chemiluminescencesignalshowedtwobrightspots,indicatingtwo ignitionpoints.Theyprovidedtwoexplanations.Thefirstone sug-geststhat the second brightspotmay bea reflectionof thefirst initial combustionsignalfromthedownstreamsurfaceofthe un-shockedbubbleinterface.Thesecondexplanationassumesthat si-multaneous detonations attwo specificpointsin thecompressed bubble are possibleas theinduction time ofthe second reaction spotisshortenedbythehighercompressionoftheshockfocusing. Apropercombinationofinductiontimesthuscanleadto simulta-neousignitions.OurnumericalresultsforMa4=2.50supportthis

latterhypothesis.

6.3. Limitationsandcriticaldiscussion

The comparison of our two-dimensional simulations with the experimental dataofHaehn etal. [28]has somelimitations with respecttothebubbleevolution.

Asshownin Fig.15,thetransverseexpansionofthebubblegas deviatesfromtheexperimentalresultsinthelong-termevolution. RecentstudiesofWangetal.[59]observethesamewithrespectto thetransversebubblediameteroftwo-andthree-dimensionalSBIs in their experiments. The vortex stretching term of the vorticity equation vanishes fora two-dimensional simulation.The missing ofthistermaffectstheexpansionandincreasesthetransverse ex-pansioninthelong-termevolution [60],whichexplainsthe devia-tionofthetransversebubblediameterfort∗6.Hejazialhosseini et al. [60] investigated thevortex dynamics inthree-dimensional inertSBIandshowedthatthegrowthrateofthevortexstretching term increases signiﬁcantly in the long-term evolution and con-tributestoadecreaseofthetransversebubblediameter,aneffect missingintwo-dimensionalsimulations.

Furthermore, the shock-focusing in three dimensions is stronger, which shortens the ignition delay time. To compensate forthis effect,we simulateRSBI ata slightly highershock Mach number to achieve the same ignition delay time. The validity of two-dimensional simulations containing shock-induced instabili-ties has also been shown by Peng et al. [61] in their study of vortex-accelerated secondary baroclinicvorticity deposition. Klein et al. [62] investigated the interaction between a sphere and a shockwaveathighshockMachnumbers.Theycomparedthe two-dimensionalresultswithexperimentaldataandobservedgood ac-cordanceintheradialandaxialwidthoftheshockedsphere.Both studies achieved very good agreement betweentwo-dimensional RMIandexperiments,eveninthelong-timedynamics.

Nevertheless some effects are not resolved by our simulation such as the onset of turbulence. Three-dimensional vortex rings tendtobecomeunstableandvortexstretchingmayeventually re-sult in broad-band turbulence [63]. This production mechanism issuppressedinatwo-dimensionalsimulation.Three-dimensional effectscannotbeneglected,howevertheybecomeimportantonly inthelong-termevolution.Niederhausetal. [64]studiedSBIsand showed that three-dimensional effects affect the total enstrophy onlyatlatetimes.Milesetal. [65]supportthisassumption,asthey alsoobserved no significantdifferences ofthe earlygrowth rates ofshock-inducedinstabilitiesbetweenthree-andtwo-dimensional simulations. Further studies report that vortex stretching affects only the long-term evolution of the mixing rate [60]. These ob-servationssupporttheintegrity ofourresults,asthechemical re-actionanditsinteractionwiththehydrodynamicinstabilitiesoccur intheearlystageofSBI.Theverygoodagreementbetweenour nu-mericalresultsandtheexperimental dataofHaehnetal. [28] in-dicate that three-dimensional effects maynot be very significant forthe specificphenomenaconsideredhere.Hence, weare confi-dent to providevalid andreliable results forthe investigation of reactionwave characteristicsandits influenceon theglobal bub-bledynamics.Inparticularthegoodagreementofthenormalized transverse bubble diameter for the detonating RSBIwith experi-mentaldata indicates that essential mechanisms are reproduced. Furthermore,thedetectionandanalysisofadoubledetonationin thesimulations supports theexperimental observationsof Haehn etal. [28].

7. Conclusion

Wehavepresentedresultsof areactingshock–bubble interac-tionatdifferentshockMachnumberswithdetailedH2–O2

chemi-cal reactionkinetics. Agas bubbleﬁlled witha reactive stoichio-metric gas mixture of H2, O2 and Xe is penetrated by a shock

wavewithMachnumbersbetweenMa₌2_.13andMa₌2_.90.The planar shock wave propagates through the domain and interacts with the cylindrical density inhomogeneity, inducing Richtmyer– Meshkovinstabilities.Theconvergentshape ofthebubblefocuses theshock,whichtriggersignitionofthebubblegas.Dependingon theshockMachnumber,thepressuresensitiveH2–O2gasmixture

showsdifferentinductiontimes,ignitionspotsandreactionwave types,which stronglyaffect the spatialbubbleevolution and the mixingprocess. Aweak shockwave induces a deﬂagration wave, highershockMach numbersdrivehigh-pressurereactions, result-inginadetonationwave.

WeshowedthatthevariationoftheshockMachnumbercovers severalreaction wave types with differentimpact on the mixing process.A deflagration wave hasa minorinfluenceon theglobal bubble evolution and leadsto a flow field dominated by hydro-dynamiceffects(Da≈ 0.28).The growthofsecondaryinstabilities ispartiallyaffected,whichdecreasesmixingbyabout30%.Higher shockMachnumbersandthesubsequentdetonationwavesleadto achemicallydrivenflowfield(Da1.32).Thesupersonicreaction