Calcolo dello Stato di Equilibrio -...

CAMERA DI COMBUSTIONE DESCRIZIONE DELLA CHIMICA Calcolo dello Stato di Equilibrio CAMERA DI COMBUSTIONE:

Calcolo dello Stato di Equilibrio


In un processo di combustione adiabatico ed isobaro in cui la combustioneraggiunge uno stato di equilibrio chimico valgono le seguenti relazioni:

Vincolo stechiometrico (ponderale): Conservazione delle moli atomiche:

Condizione termodinamica di equilibrio chimico(max entropia == min entalpia libera di Gibbs)

Conservazione dell’energia in termini di entalpia assoluta della miscela

Saturday, July 11, 15

TURBOGETTO SEMPLICE CICLO A PUNTO FISSO CICLO IN VOLO Ese. 5: Turbo-Getto TURBOFAN TURBOF

COMBUSTORE

ANALISI A PUNTO FISSO: COMBUSTORE

EQ. ENERGIA (Ls = 0, M ≪ 1): ∆h0 ≃ ∆h = ∆Q

ηb := Q/(mf Qf ) RENDIMENTO DI COMBUSTIONE

T4: TEMPERATURA ALL’INGRESSO DELLA TURBINA (TIT)

ηpb := p4/p3 RENDIMENTO PNEUMATICO DEL COMBUSTORE

f := mf/ma RAPPORTO COMBUSTIBILE/ARIA, O DI DILUIZIONE

mah3 + mfhf + Q = (ma + mf )h4

f ≪ 1 ⇒ mah3 + Q = mah4 ⇒ macp(T4 − T3) = Q

Q = ηb mf Qf

f =cp (T4 − T3)

ηb Qf=

cpT3

ηb Qf

!

T4

T3− 1

"

p4 = ηpb p3



Processi chimicamente reversibili: congelati, in equilibrio

Processi chimicamente reversibili: congelati, in equilibrio

Irreversible processes

dS = dSext + dSint

dSext =1

TdU +

p

TdV =

dQext

T= 0 Non adiabatic system

dSint = −1

T

!

j

µjdNj > 0 Chemically reactive system

Internal (chemical) reversible processes

dSint = −1

T

!

j

µjdNj = 0 ⇔ dG =!

i

µidNi = 0

Chemically frozen processes (air intake, compressor, turbine, nozzle)

∀j : dNj = 0 ⇒ Nj = const ⇒ c(v,p)(T, Yi) = const if gas is calorically perfect

Processes in chemical equilibrium (combustion chamber)

dG =!

j

µjdNj = 0 ⇒ Minp,Tgiven

[G (p, T,Nj)] ⇒ Nj = N∗

j (p, T )



Conservazione delle moli atomiche: stechiometria


Supponiamo che la miscela sia formata da 8 specie chimiche:

1 2 3 4 5 6 7 8H2 O2 H O OH H2O HO2 H2O2

La conservazione delle moli atomiche si esprime con 2 equazioni algebriche

NH = 2n1 + n3 + n5 + 2n6 + n7 + 2n8

NO = 2n2 + n4 + n5 + n6 + 2n7 + 2n8

che in forma matriciale si scrive

Adn =

(2 0 1 0 1 2 1 20 2 0 1 1 1 2 2

)

⎧

⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨

⎪⎪⎪⎪⎪⎪⎪⎪⎪⎩

dn1

dn2

dn3

dn4

dn5

dn6

dn7

dn8

⎫

⎪⎪⎪⎪⎪⎪⎪⎪⎪⎬

⎪⎪⎪⎪⎪⎪⎪⎪⎪⎭

= 0

questo sistema lineare di 2 equazioni in 8 incognite ammette ∞6 = 8 − 2 soluzioni.





Le ∞6 soluzioni si trovano partizionando la matrice A ed il vettore dn:

(2 00 2

){dn1

dn2

}

= −

(1 0 1 2 1 20 1 1 1 2 2

)

⎧

⎪⎪⎪⎪⎪⎨

⎪⎪⎪⎪⎪⎩

dn3

dn4

dn5

dn6

dn7

dn8

⎫

⎪⎪⎪⎪⎪⎬

⎪⎪⎪⎪⎪⎭

per ogni scelta delle componenti dnj , j = 3, 8 si trovera una ed una sola coppiadnj , j = 1, 2 che soddisfa il sistema di due equazioni

come si possono scegliere le dnj , j = 3, 8 in modo di essere sicuri di averle presetutte ? si utilizza una base di vettori linearmente indipendenti ovvero:

ek ={

eik

}

=

{1 se i = k0 se i = k

}

k = 1,8 - 2 = 6

Con la scelta degli ek effettuata utilizzando vettori linearmente indipendenti

consente di scrivere le ∞6 soluzioni in questo modo:

⎧

⎪⎨

⎪⎩

{dn1

dn2

}

= −

(2 00 2

)−1 (

1 0 1 2 1 20 1 1 1 2 2

)

dn2+k

dn2+k = ekdξk

k = 1, 6





Tornando ad una rappresentazione per componenti si ottiene

dn1 = − 12dξ

1

dn2 = 0dn3 = dξ1

⎫

⎬

⎭⇒

1

2H2 = H dξ

1assume valori compresi ≥ 0

che indica che ogni ∞1 soluzioni rappresenta una reazione chimica virtuale

E qundi in una miscela con 8 specie formata da 2 elementi atomici si possono avereal piu 8-2=6 reazioni chimiche linearmente indipendenti

in forma compatta si puo scrivere

dni =∑

k

νki dξ

k

Le 6 direzioni cosı trovate individuano un sottospazio in R8 che in algebra lineare

viene chiamato: spazio nullo della matrice A ( NullSpace[A] in Mathematica ).

Lo spazio nullo e il sottospazio in cui, a partire da una composizione iniziale dellamiscela assegnata, le reazioni chimiche trasformano la miscela in modo che ilnumero di atomi iniziali si conservi

dn1dt = − 1

2dξ1

dt

dn2dt = 0

dn3dt = dξ1

dt

dξ1

dt:= r

1(p, T,Nj) = r1f − r

1b



Condizione termodinamica di equilibrio chimico

Condizione termodinamica di equilibrio chimico

Reversible (zero entropy) process dG =∑

j

µjdNj = 0

Stoichiometric constraint dNj =∑

k

νkj dξk

⇓

dG =∑

j

µj

∑

k

νkj dξk =∑

k

dξk∑

j

µjνkj = 0

⇓

∀dξk = 0,

⎧

⎨

⎩

∑

j

µjνkj = 0

⎫

⎬

⎭

Nreactions

k=1

⇓

Free Enthalpy (Gibbs) is stationary ⇔ Chemical equilibrium



Legge di Azione di Massa

Legge di Azione di Massa

Equilibrium Condition∑

j

µj (T, pj) νkj = 0 k = 1,Nreactions

µj(T, pj) := Hj(T ) − T Sj(T, pi) Sj(T, pj) := S0j (T ) − ℜLog(

pj

pref

)

µj(T, pj) := Hj(T ) − T

[

S0j (T ) − ℜLog(

pj

pref

)

]

=(

Hj(T ) − T S0j (T )

)

+ ℜTLog(pj

pref

) = µ0j (T )

∑

j

[

µ0j (T ) + ℜTLog(

pj

pref

)

]

νkj = 0 ⇒

∑

j

[

µ0j (T )

]

νkj = −

∑

j

[

ℜTLog(pj

pref

)

]

νkj

∑

j

[

Hj(T ) − T S0j (T )

]

νkj = −ℜT

∑

j

νkj Log(

pj

pref

) = −ℜT∑

j

Log(pj

pref

)νkj

Exp

⎧

⎨

⎩−

1

ℜT

∑

j

[

Hj(T ) − T S0j (T )

]

νkj

⎫

⎬

⎭

︸︷︷︸

Kp(T )

= Exp

⎧

⎨

⎩

∑

j

Log(pj

pref

)νkj

⎫

⎬

⎭

Kkp (T ) = Exp

⎧

⎨

⎩Log

∏

j

(pj

pref

)νkj

⎫

⎬

⎭=

∏

j

(pj

pref

)νkj

Law of Mass Action Kkp (T ) =

∏

j

(pj

pref

)νkj





Ns + 2 incognite : composizione della miscela (N jProducts, j = 1, Ns), Numero

totale di moli (Ntot), e temperatura adiabatica di fiamma (TProducts)

Ns + 2 equazioni:Conservazione delle moli atomiche:

Ni =

Ns!

j=1

aijN

jReac =

Ns!

j=1

aijN

jProducts i = Ne

Condizione termodinamica di equilibrio chimico (Legge di azione di massa)

Kkp (T ) =

"

j

(pj

pref

)νkj k = Ns − Ne pj =

Nj

Ntot

p

Numero totale di moli

Ntot =Ns!

j=1

Nj

Conservazione dell’energia in termini di entalpia assoluta della miscela

H(TReac, NjReac) = H(TProducts, N

jProducts)



Risultati del Calcolo dello Stato di Equilibrio

Risultati del Calcolo dello Stato di Equilibrio

MISCELA H2/O2

0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0

2900

3000

3100

3200

3300

3400

Figure: Variazione Temperatura adiabaticadi equilibrio con rapporto di equivalenza(p=1 e 10 atm); T reagenti = 300K.

0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0

0.001

0.005

0.010

0.050

0.100

0.500

1.000

Figure: Variazione Composizione diequilibrio con rapporto di equivalenza (p=1e 10 atm); T reagenti = 300K.



FLUSSO ALLA RAYLEIGH

FLUSSO ALLA RAYLEIGH

Variazione di T0 e variazione di Mach in CC

!

T0

T⋆0

"

cal

[M, γ] =2M2

#

1 + δM2$

(γ + 1)#

1 + γM2$2

!

T04

T03

"

cal

[M3,M4, γ] =

%

T0T⋆

0

&

cal[M4, γ]

%

T0T⋆

0

&

cal[M3, γ]

=

#

1 + δM24

$ #

M4 +M23M4γ

$2

#

1 + δM23

$ #

M3 +M3M24γ

$2

Variazione di p0 e variazione di Mach in CC

!

p0p⋆0

"

cal

[M, γ] =2#

1 + δM2$

γγ−1 (γ + 1)

11−γ

1 + γM2

!

p04

p03

"

cal

[M3,M4, γ] =

%

p0p⋆0

&

cal[M4, γ]

%

p0p⋆0

&

cal[M3, γ]

=

!

1 + δM23

1 + δM24

"

γγ−1 1 + γM2

4

1 + γM23



RENDIMENTO PNEUMATICO

RENDIMENTO PNEUMATICO

Rapporto T04/T03 in funzione della portata di combustibile f

f =cp (T4 − T3)

ηb Qf

≈cpT03

ηb Qf

!

T04

T03− 1

"

⇒T04

T03= 1 +

f ηb Qf

cpT03

Numero di Mach in ingresso M3 della CC per assegnati Mach in uscita M4 erapporto T04/T03

Solve

#

T04

T03==

$

2 +M24 (γ − 1)

% $

M4 +M23M4γ

%2

$

2 +M23(γ − 1)

% $

M3 +M3M24 γ

%2 ,M3

&

Perdite di pressione p0 in funzione della variazione di Mach tra ingresso eduscita della camera

ηpb :=

!

p04p03

"

cal

[M3,M4, γ] =

!

1 + δM23

1 + δM24

"

γγ−1 1 + γM2

4

1 + γM23


FP7-SPACE-CALL-1

In Space PropulsionISP-1

Space Propulsion Conference – Bordeaux – May 7-10, 2012

Micro Combustion ChamberM3 Test Bench - DLR

13

M3 chamber characterizes ignition behaviour, with varying injection geometries and flow conditions

Mass flows are determined through Coriolis flow meters

Sonic nozzles set the correct mass flows

Sonic nozzles provide an effective separation between feed lines and injector head, thus minimising low frequency combustion instabilities


ISP-1 - September 22-23, 2011

Model Equations for Unsteady CSTR with Energy Deposition

d(ρV )dt

= !mfuel[t,p]+ !moxid[t,p]− !mnozzle[p,T,Y]

d(ρhV )dt

= !mfuel[t,p]hfuel + !moxid[t,p]hoxid − !mnozzle[p,T,Y]h[T,Y]+Vdpdt

+ !q[t]V

dYj

dt=

1ρ[p,T,Y]V

!mfuel[t,p]Yj + !moxid[t,p]Yj − !mnozzle[p,T,Y]Yj( ) + Wj

ρ[p,T,Y]ω j[p,T,Y]−

Yj

ρ[p,T,Y]dρdt[p,T,Y]

p = ρRT →1ρdρdt

=1pdpdt

−1RdRdt

−1TdTdt

R =RUWjj

∑ Yj →dRdt

=RUWjj

∑ dYj

dt

h = hjj∑ Y →

dhdt

= hjj∑ dYj

dt+ Cp

dTdt

→ Cp = Cj , pYjj∑

!q(t)[ ErgCentimeter3Second

] = ε[Joule]4πσ r

3σ t

e−12t− t0σ t210[ Erg

Joule]

!mfuel[t,p] = !mred[p

p0,fuel (t), ηfuel , γ fuel ]

p0,fuel (t)Afuel

RfuelT0,fuel

!moxid[t,p] = !mred[p

p0,oxid (t), ηoxid , γ oxid ]

p0,oxid (t)Aoxid RoxidT0,oxid

!mnozzle[p,T,Y] = !mred[pa

p, ηoxid , γ (T ,Y )] p Aoxid

R(Y)T

!mred[Πp , η, γ ] = 2γγ −1

1Πp

η 1− Πp

γγ −1

⎛

⎝⎜

⎞

⎠⎟

1−η 1− Πp

γγ −1

⎛

⎝⎜

⎞

⎠⎟

with Πp ≥2

γ +1⎛⎝⎜

⎞⎠⎟

γγ −1



Model Equations for Unsteady CSTR with Energy Deposition

dYj

dt=

R[Y ]TpV

!mfuel[t,p]Yj + !moxid[t,p]Yj − !mnozzleYj( )

+Wjω j[p,T ,Y ]ρ[p,T ,Y ]

−Yj

ρ[p,T ,Y ]dρdt

[p,T ,Y ] j = 1 , Ns

dTdt

=- !mnozzle[p,T ,Y ] + !mfuel[t,p]

hfuel − h[T ,Y ]R[Y ]T

+ 1⎛⎝⎜

⎞⎠⎟

+ !moxid[t,p] hoxid − h[T ,Y ]R[Y ]T

+ 1⎛⎝⎜

⎞⎠⎟

(Cp-R) p V

−hj[T ]

j∑ dYj

dt -T RU

Wj

dYj

dtj∑

(Cp[T ,Y ]-R[Y ])+

!q[t]ρ[p,T ,Y ] (Cp[T ,Y ]-R[Y ])

dpdt

=p γ [T ,Y ]R[Y ]T

- γ [T ,Y ]−1γ [T ,Y ]

hj[T ]j∑ dYj

dt + T RU

Wj

dYj

dtj∑

⎛

⎝⎜⎞

⎠⎟

+ γ [T ,Y ]R[Y ]TV

- !mnozzle[p,T ,Y ] + !mfuel[t,p]hfuel − h[T ,Y ]Cp[T ,Y ] T

+ 1⎛⎝⎜

⎞⎠⎟

+ !moxid[t,p] hoxid − h[T ,Y ]Cp[T ,Y ] T

+ 1⎛⎝⎜

⎞⎠⎟

⎛⎝⎜

⎞⎠⎟

!q[t][ ErgCentimeter3Second

] = ε[Joule]4πσ r

3σ t

e−

12t− t0σ t

210[ Joule

Erg]

R[Y ] = RUWj

Yjj∑

ρ[p,T ,Y ] = pR[Y ]T

h[T ,Y ] = hj[T ]j∑ Yj

dρdt

=1V!mfuel[t, p]+ !moxid[t, p]− !mnozzle[p,T ,Y ]( )

!mfuel[t,p] = !mred[p

p0,fuel (t), ηfuel , γ fuel ]

p0,fuel[t]Afuel

RfuelT0,fuel

!moxid[t,p] = !mred[p

p0,oxid (t), ηoxid , γ oxid ]

p0,oxid[t]Aoxid RoxidT0,oxid

!mnozzle[p,T,Y] = !mred[pa

p, ηoxid , γ (T ,Y )] p Aoxid

R[Y]T

!mred[Πp , η, γ ] = 2γγ −1

1Πp

η 1− Πp

γγ −1

⎛

⎝⎜

⎞

⎠⎟

1−η 1− Πp

γγ −1

⎛

⎝⎜

⎞

⎠⎟

with Πp ≥2

γ +1⎛⎝⎜

⎞⎠⎟

γγ −1



Data

p0fuel = 2.5 patm; pfuel = p0fuel;Tfuel = 290 Kelvin; T0fuel = Tfuel;ηfuel = 0.8;

p0oxid = 5.0 patm; poxid = p0oxid;Toxid = 290 Kelvin; T0oxid = Toxid;ηoxid = 0.8;

radiusoxid = 0.11 Centimeter; Aoxid = Pi radiusoxid^2;innerradiusfuel = 0.24 Centimeter;outerradiusfuel = 0.30 Centimeter; Afuel = Pi (outerradiusfuel^2 - innerradiusfuel^2);

Lchamber = 16 Centimeter;radiuschamber = 3 Centimeter;Volume = Lchamber radiuschamber;

pa = 1 patm; p0a = pa;Ta = 290 Kelvin; T0a = Ta;radiusnozzle = 1 Centimeter;Anozzle = Pi radiusnozzle^2;ηnozzle = 0.8;

T0 = 300. Kelvin; p0 = 1. patm;



Inputs

0.0 0.5 1.0 1.51¥106

2¥106

3¥106

4¥106

5¥106

0.0 0.5 1.0 1.50

5.0¥108

1.0¥109

1.5¥109

2.0¥109

2.5¥109

3.0¥109

Total Pressures at Fuel & Oxygen Manifold in time

Energy deposition in time

sigr = 1.0*10^-3 Meter;sigt = 2.5*10^-4 Second;eps = minimumEnergy ;c1 = eps/(4 Pi sigr*sigr*sigr*sigt);qdotSI = c1*Exp[-.5*((time - t0)/sigt)^2];SItoCGS = 10;qdot = qdotSI * SItoCGS



Outputs: Methane/Oxygen SystemSensitivity of peak pressure on energy deposition time

Temperature in time

Peak pressure vs time lag

Pressure in time

• short energy deposition time allow little reactants to fill the chamber -> lower total reactant mass yields a lower pressure peak during ignition

• short/long energy deposition time: ignition develop three/four stages

• first stage[heating]: isobaric T growth due to energy addition• second stage[explosion] (kinetics >> convection): T increases due to fast kinetics, and P follows for inertial confinement

because the kinetics is faster than convection• third stage[relaxation] (kinetics ~ convection): T&P decrease because the higher P produces a larger outflow• fourth stage[near-equilibrium] (kinetics << convection): T&P increase following the arrival of reactants which are instantly

burned

0.10 0.15 0.20 0.25 0.305001000150020002500

0.10 0.15 0.20 0.25 0.301.0¥1061.2¥1061.4¥1061.6¥1061.8¥106

0.5 1.0 1.51.0¥106

1.2¥106

1.4¥106

1.6¥106

1.8¥106



Outputs: Methane/Oxygen SystemSensitivity of peak pressure on minimum energy deposition

Temperature in timePeak pressure vs minimum energy deposition

Pressure in time

• too low power levels cannot ignite the mixtures • higher power levels shorten the ignition times• less reactants can fill the chamber• peak pressure exhibits a max at ~50mJ

0.15 0.16 0.17 0.18 0.19 0.20

500

1000

1500

2000

2500

3000

3500

0.15 0.16 0.17 0.18 0.19 0.20

1.1¥106

1.2¥106

1.3¥106

1.4¥106

1.5¥106

1.6¥106

0.02 0.04 0.06 0.08 0.10 0.12 0.141.58¥106

1.59¥106

1.60¥106

1.61¥106

1.62¥106

1.63¥106


FP7-SPACE-CALL-1



To check the ability of the physical models and prediction tools to reproduce:

• ignition model • flame propagation

- kernel formation- flame kernel evolution in a turbulent flow

• anchoring process

A dedicated experiment has been carried out at the M3 Test Bench (DLR Lampoldhausen) and has been used as a reference test case to benchmark different numerical approaches

Ignition is triggered using a laser beam to control the ignition point location and energy release, in a well controlled gas/gas injection configuration

This way a clean experimental configuration is obtained, allowing to check the ability of the numerical tools to reproduce flame propagation and anchoring

Two test campaigns have been carried out:

• ambient pressure• low pressure

Objectives of ISP-1 WP 2.5

20


FP7-SPACE-CALL-1


Space Propulsion Conference – Bordeaux – May 7-10, 2012 21

Test case computed : 21.07.2011-7 (ambient, attached)

Ignition Sequence

- 117 -57 0 ignition Time line ms

The ignition overpressure (if all the quantity of propellant in the chamber prior to ignition is burnt instantaneously) can be found from:

pmaxpc

=!miτ i!mτ r

Decreasing the ignition delay decreases the overpressure during ignition.

τ r :=

M!me

These findings helped the selection of an ignition sequence of the test campaign able to minimize the overpressure.

An accurate selection of the ignition sequence is essential to avoid undesired peaks in the chamber pressure.

τ i is the ignition delay

τ r is the mean residence time in the chamber

!me is the massflow rate of gases exhausted though the nozzle

M is the mass of gas inside the chamber

!miτ i is the mass of reactants in the chamber over the time τ i

!mτ r is the mass of propellants in the chamber at any time

Accumulated unburnt masses prior to ignition determine the pressure rise during ignition and its overpressure.


FP7-SPACE-CALL-1



Gaseous Oxygen, before ignition

Coaxial Methane, before ignition, with Oxygen

Schlieren images of laser pulse (150 mJ); Expanding blast wave; Sequence sampled from multiple laser shots therefore only an estimate for the timing can be given; Time step ~ 50 µs

Schlieren image sequence of laser ignition of CH4/O2 gas/gas

COLD FLOW

HOT FLOWLASER-PULSE ON QUISCENT NITROGEN

Micro Combustion ChamberM3 Experiment - DLR

22


FP7-SPACE-CALL-1



Main experimental findings from DLR

23

Flame lift-off length inversely proportional to chamber pressure

Excitation caused by heat release results in spontaneous emission of OH radicals at ca. 305-309 nm

Spontaneous emission in the UV-range is recorded by an intensified high-speed CCD video camera

Complete chamber is visualised with a resolution of 512 x 256 pixels

A band-pass filter (310 nm ± 5 nm) selects only OH-radical emission

Ambient pressure ignition

CC pressure rising due to heating


FP7-SPACE-CALL-1



Blast Wave

Motivations for CFD Analyses

Recirculation

Fuel Blockage

Quantitative analysis of ignition can be done by a "Well-Stirred Reactor" (WSR) or "Continuously Stirred Tank Reactor" (CSTR) models

Both models assume infinitely fast and efficient mixing in the chamber

These models allow to readily carry out all the relevant sensitivity analyses

Why then making the costly and tedious CFD analyses ?

Because of the critical role of multi-dimensional phenomena !!


FP7-SPACE-CALL-1



URANS MODELLING OPTIONS (CFD++ METACOMP Tech)• Flow geometry 2D Axi-symmetric • real gas, compressible equations• Ideal gas equation of state• multi-species with frozen/active detailed chemical kinetics• viscous flow• transient integration (time accurate, point implicit)• second order space discretization• turbulent modelling on: RANS two-eqns k-epsilon • Ignition by « hot spot »

ICs: • quiescent Nitrogen fills the chamber at T=290 K, and p=101325 Pa

BCs:• Walls are assumed isothermal, post tip and nozzle are considered

adiabatic• Inflow: constant total pressure and temperature for both fuel and

oxygen• Outflow: subsonic flow with prescribed ambient pressure

MESH: ~ 230K cells; block-structured

CHEMICAL KINETIC MECHANISM FOR CH4/O2+INERT NITROGEN: • GRI3.0 (53 species) with Nitrogen kinetics removed involves 36

species• Mechanism simplification (in-house tools) trimmed the mechanism to

15 species

URANS & LES Modelling Options

25

LES MODELLING OPTIONS (CEDRE ONERA)• Flow geometry 3D • real gas, compressible equations• Ideal gas equation of state• multi-species with a frozen/active global reaction mechanism• viscous flow• transient integration (time accurate)• second order space discretization• Smagorinski model for LES subgrid closure• Ignition by « hot spot »

ICs: • quiescent Nitrogen fills the chamber at T=290 K, and p=101325 Pa

BCs:• Walls are assumed adiabatic• Inflow: ramped total pressure and temperature for both fuel and

oxygen• Outflow: subsonic flow with prescribed ambient pressure

MESH: ~ 10M cells; unstructured

CHEMICAL KINETIC MECHANISM FOR CH4/O2+INERT NITROGEN: • Global kinetics (Jones and Lindstedt, adapted by Kim)


FP7-SPACE-CALL-1


Space Propulsion Conference – Bordeaux – May 7-10, 2012 26

Test case computed : 21.07.2011-7 (ambient, attached)• Step 1 : O2 cold flow– Not the whole injector is meshed : to account for pressure losses in Ox tube, boundary condition pressure is Pi=9 bars– To account for pressure increase vs time, relaxation is activated on the boundary condition– From t=-117ms to -57ms

• Step 2 : CH4 cold flow– Along with stabilized O2 cold flow– From t=-57ms to t=0ms

CFD AnalysesInflow BCs for LES

The filling process of the chamber by the cold reactants needs to be accurately replicated because the total mass of reactants at time T0 sets the total level of energy available for combustion during the ignition start-up


FP7-SPACE-CALL-1



How we trim the kinetics

Detailed (NS36) vs Simplified (NS15)

27

GRI-Mech 3.0 (53 spcs and 325 reversible rcns) is used as reference mechanism

All N-containing species are removed, except N2, together with all N-related reactions, to yield a detailed mechanism with 36 spcs and 219 reversible rcns

Mechanism simplification done by in-house tools

A spatially homogeneous, iso-choric, adiabatic, forced ignition with gaussian energy deposition drives the simplification procedure

The selected simplified mechanism involves 15 species and 57 reactions

The fastest time scale of the simplified mechanism is two orders of magnitude larger than the one of the detailed (stiffnes reduction)


O⇣

1pdpdt

⌘⇡ O

�1T

dTdt

�+O

�1R

dRdt

�+O

⇣mF

⇢V + mOx

⇢V � mOut

⇢V

⌘

FP7-SPACE-CALL-1



Test Case#07 Ambient Pressure

Chamber pressure time evolution

URANS & LES return similar chamber pressure histories

Both feature a larger-than-experiment chamber pressure peak value and growth rate

Chamber pressure peak exceeds fuel manifold pressure -> fuel blockage

Late pressure evolution exhibits oscillations tuned with chamber acoustics

NB: Chamber pressure growth rate proceeds as:


FP7-SPACE-CALL-1



Test Case ValuesMethane = 0.60 g/sOxygen = 2.26 g/s


Massflows

Methane and oxygen inflows not fixed at their choked values when the pressure drop falls below the critical value

Hot products outflow eventually fixed by sonic condition at nozzle

Nearly constant massflows at late times suggest a nearly steady condition

URANS & LES fluxes are quite in agreement one to another (despite the slight different treatment of the inflow BCs)


FP7-SPACE-CALL-1




Propellants Mach Number

Oxygen stream choked during cold flow injection (M>1)Oxygen stream subsonic when chamber pressure peaks (M<1)Methane stream never choked (M<1); tends to zero near chamber pressure peak


FP7-SPACE-CALL-1



• A blast wave propagates spherically outward, reflects at the injector plate and at the wall, head-on collides at the symmetry axis, propagates downstream towards the nozzle only to be reflected backward towards the injector plate

• The under-expanded oxygen jet gradually fades away when the flow becomes subsonic• Note the formation of transient traveling pressure peaks at the symmetry axis


Blust Wave following Laser Pulse


FP7-SPACE-CALL-1



HCO

Temperature


Kernel initiation and flame propagation

• The flame kernel propagates downstream near the symmetry axis as convected by the fast central jet

• The flame kernel propagates across the recirculation region away from the axis

• When methane is not entering the chamber anymore, the cold oxygen jet is not consumed and leaves the chamber unburned-> the temperature field at the axis becomes very cold

• A significant amount of hot products is still present in the chamber ready to reignite the propellants when the chamber pressure is lowered by the mass loss through the choked nozzle


FP7-SPACE-CALL-1



• Methane blockage is clearly monitored by this movie• Note the role of the cavity as an accumulator of the blocked methane• The blockage process exhibits fluctuations coupled with the chamber acoustics

Methane

Oxygen


Fuel blockage


FP7-SPACE-CALL-1



LESTemperature iso-contour at T=2000K

Overall LES dynamics consistent with URANS predictions

LES captures the non symmetric flow evolution triggered by the off-axis laser pulse location

LES captures 3-D jet instabilities

LES captures a “breathing” evolution of the flame front as also noticed in experiments

LES captures the development of cold unreacted pockets


FP7-SPACE-CALL-1



URANS OH vs Experimental OH*

The comparison shows that the computed OH field is qualitatively exhibiting a similar shape and shape evolution, moves downstream at about the right speed, and posses a brighter core region surrounded by a darker halo, where locally the light intensity is proportional to the amount of OH


FP7-SPACE-CALL-1



LES Temperature vs Experimental OH*

The shape and position of the OH∗ emission field is quite similar to the one of the temperature and reactive zones in the LES computation

OH* emission field accounts for the integral contribution along the chamber cross-section, whereas both LES and URANS results refer to cut-planes passing through the chamber axis


FP7-SPACE-CALL-1



Flame Index at “stationary” conditions

FI = ∇yFuel ⋅∇yOxidFI > 0 ⇒ fuel and oxidizer fluxes aim at same direction (pre-mixed nature) FI < 0 ⇒ fuel and oxidizer fluxes aim at opposite directions (non pre-mixed)

(nearly) stationary flame is: lifted, and pre-mixedCo-axial Jet Mixing layer

Pre-heat, pre-mixed flame region(Solid lines are HCO mass fraction lines)

Fuel and oxidizer from both jet and recirculation region mix here

Negative T IndexNon-premixed mixture

Positive T IndexPremixed mixture

llift−off ∼ vjetτ ignvjet ∼ p0 − pCCllift−off ∼ p0 − pCCτ ign


FP7-SPACE-CALL-1



Main problem found• Larger-than-experiment chamber pressure peak value and growth rate

Possible causes (given that the pressure peak is mostly linked to the mass of CH4 and O2 in the chamber at T0)

• neglecting nitrogen filling the propellant manifold at T0 (volume of pipe(s) between the probe and the boundary conditions) realizes "too much" reactants in the chamber

• supersonic oxygen jet spreads too quickly by numerical dissipation and causes an excess of oxidant in the chamber at T0

Conclusions

Lesson learned

• URANS axi-symmetric calculations can be effectively able to provide a rather detailed picture of the ignition events, albeit there remains a number of issues for the quantitative accuracy of the URANS predictions

• 2D axi-symmetric URANS and 3-D LES provide predictions in satisfactory agreement, even when rather different kinetic mechanisms have been adopted

• CFD analyses offered interesting contributions in the understanding of a number of critical ignition phenomena, which are difficult to appreciate on the basis of experimental diagnostics alone


FP7-SPACE-CALL-1



This work has been carried out with the support of:

FP7 EU Grant no.218849, titled "In-Space Propulsion-1" (ISP-1)

M.Valorani acknowledges the support of:

CASPUR Competitive HPC Grant 2009

The URANS flow solver is: CFD++ by Metacomp Technologies, Inc.

The LES flow solver is: CEDRE, an ONERA in-house software package

Acknowledgements

39


Calcolo dello Stato di Equilibrio -...

Documents

Transcript of Calcolo dello Stato di Equilibrio -...