TEXT(CO & H2 Diffusion Flame
TITLE
DISPLAY
The case considered is the turbulent combustion of CO and H2 with
O2 in a pipe of 32cm diameter. Concentric 'fuel' and 'oxidiser'
streams enter the pipe axially, and then having mixed and burnt,
exit axially 2m downstream of the inlet. The fuel stream, which
comprises H2, CO, H2O, CO2 and N2, enters through an inner
injector of 8cm diameter, and the oxidiser stream, which comprises
O2, N2, H2O and CO2 enters through an outer annulus. The tempera-
tures of the fuel and oxidiser streams are 1754 and 1518 Deg C,
and the corresponding mass flow rates are 0.0726 and 0.144 kg/s.
The fast-chemistry assumption is invoked for main combustion
calculation in terms of two competing single-step infinite-rate
reactions: 2H2 + O2 > 2H2O and CO + O2 > CO2. The NOX calculation
is peformed as a post-processing run with a frozen hydrodynamic
and temperature fields. The NO distributions are calculated by
CHEMKIN via a reaction mechanism involving the Zeldovich reactions
and fuel kinetics.
ENDDIS
There are 2 runs to be done:
Run 1 - use the generalised SCRS to do the main combustion
problem (RUN1=T) plus CHEMKIN to calculate the
temperature from the enthalpy ( this requires the
file SCRS.CKM for CHEMKIN ).
Run 2 - restart with frozen flow, temperature and turbulence
fields and solve the detailed chemistry with NOx
included using the implicit solver in the CHEMKIN
Interface (RUN1=F). This requires the file RADN.CKM
for CHEMKIN.
The inlet conditions for the case are:
Fuel Stream Oxidiser Stream
M = 0.0726 kg/s M = 0.144 kg/s
W = 96.0 m/s W = 9.9 m/s
T = 1754 C T = 1518 C
YCO = 0.091 mol/mol YO2 = 0.086 mol/mol
YH2 = 0.089 " YCO2= 0.055 "
YCO2= 0.036 " YH2O= 0.111 "
YH2O= 0.164 " YN2 = 0.748 "
YN2 = 0.62 "
The expected outlet mass fractions for the SCRS run are:
YO2 = 0.026 YCO2 = 0.133 YH2O = 0.732 YN2 = 0.109
GROUP 1. Run title and other preliminaries
TEXT(RADIAN: 2D Turbulent Diffusion Flame
CHAR(CPDF);BOOLEAN(RUN1,HSOLV,THRAD,BLOCK)
HSOLV=T;THRAD=T
MESG( porosity test ? (default=N)
READVDU(ANS,CHAR,N)
IF(:ANS:.EQ.Y) THEN
+ BLOCK=T
+ MESG( Porosities present
ELSE
+ BLOCK=F
+ MESG( No porosities
ENDIF
MESG( Do you want a NOX run? (default=N)
READVDU(ANS,CHAR,N)
IF(:ANS:.EQ.Y) THEN
+ RUN1=F
+ MESG( NOX run selected.
ELSE
+ RUN1=T
+ MESG( Main combustion run selected.
ENDIF
REAL(WINF,WINO,KEINO,EPINO,KEINF,EPINF,PRADO,PRADI,CD)
REAL(TEMFU,TEMOX,TWAL)
TWAL=100.+273.
PRADO=0.16;PRADI=0.04;CD=0.1643;WINF=96.0;WINO=9.9
KEINO=(0.05*WINO)**2; EPINO=CD*KEINO**1.5/(0.1*(PRADO-PRADI))
KEINF=(0.05*WINF)**2; EPINF=CD*KEINF**1.5/(0.1*PRADI)
TEMFU=1754.+273.;TEMOX=1518.+273.
REAL(YFU1IN,YOX1IN,YN21IN,YH21I,YCO1I,YH2O1I,YCO21I)
REAL(YFU2IN,YOX2IN,YN22IN,YH22I,YCO2I,YH2O2I,YCO22I)
** fuel stream
YFU1IN=0.0;YOX1IN=0.0;YH21I=-0.089;YCO1I=0.091;YH2O1I=0.164
YCO21I=0.036;YN21IN=1.-YFU1IN-YOX1IN+YH21I-YCO1I-YH2O1I-YCO21I
** oxidiser stream
YFU2IN=0.0;YOX2IN=0.086;YH22I=0.0;YCO2I=0.0;YH2O2I=0.111;
YCO22I=0.055;YN22IN=1.-YFU2IN-YOX2IN-YH22I-YCO2I-YH2O2I-YCO22I
GROUP 3. X-direction grid specification
CARTES=F;XULAST=0.1
GROUP 4. Y-direction grid specification
INTEGER(NYF,NYO,NYG);NYF=4;NYO=4;NYG=NYF+NYO
IF(BLOCK) THEN
+ NREGY=3;NY=12
ELSE
+ NREGY=2;NY=8
ENDIF
IREGY=1;GRDPWR(Y,NYF,0.04,1.0);IREGY=2;GRDPWR(Y,NYO,0.12,1.0)
IF(BLOCK) THEN
+ IREGY=3;GRDPWR(Y,4,0.04,1.0)
ENDIF
GROUP 5. Z-direction grid specification
NZ=8;GRDPWR(Z,NZ,2.0,1.5)
MESG( BFC test ? (default=N)
READVDU(ANS,CHAR,N)
IF(:ANS:.EQ.Y) THEN
+ BFC=T;NONORT=T
+ MESG( BFCs present
ELSE
+ BFC=F
+ MESG( No BFCs
ENDIF
GROUP 7. Variables stored, solved & named
IF(RUN1) THEN
+ SOLVE(P1,V1,W1)
ELSE
+ STORE(P1,V1,W1)
+ NAME(15)=ENT2;NAME(45)=F;STORE(F)
ENDIF
STORE(VIST,DEN1,TMP1);SOLUTN(P1,P,P,Y,P,P,P)
SOLUTN(V1,P,P,P,P,P,N);SOLUTN(W1,P,P,P,P,P,N)
IF(HSOLV) THEN
+ SOLVE(H1);TERMS(H1,N,P,P,P,P,P)
+ SOLUTN(H1,P,P,Y,P,P,P)
ELSE
+ STORE(H1)
ENDIF
** store specific heat
STORE(SPH1,YSUM)
** activate radiation
IF(THRAD) THEN
+ REAL(ABSORB,SCAT,SIGMA,EMPW,EMISW,EMISG,EMPG)
+ ABSORB=0.35;SCAT=0.0; EMISG=0.35; EMISW=1.0
+ SIGMA=5.6697E-8; EMPW=SIGMA*TWAL**4; EMPG=SIGMA*EMISG
+ RADIAT(RADI,ABSORB,SCAT,H1)
+ SOLUTN(SRAD,P,P,Y,P,P,P);SOLUTN(H1,P,P,Y,P,P,P)
ENDIF
IF(RUN1) THEN
*** Start of extended SCRS model settings, which are equivalent to,
and supersede, the SPEDAT statements which have been de-activated
+ PRESS0=1.0000E+05
+ INTEGER(NSPEC,NELEM);NSPEC=7;NELEM=4
+ INTEGER(NCSTEP,NCREAC);NCSTEP=1;NCREAC=2
+ SCRS(SYSTEM,NCSTEP,NCREAC,NELEM,FASTC)
+ SCRS(SPECIES,CH4,O2,H2,CO,H2O,CO2,N2)
+ SPEDAT(SET,SCRS,NRSYS,I,-3)
** Define fuel, oxidiser, & intermediates
+ SPEDAT(SET,SCRS,CFUEL,C,CSG5);SPEDAT(SET,SCRS,COXID,C,O2 )
+ SPEDAT(SET,SCRS,CFP1 ,C,CO );SPEDAT(SET,SCRS,CFP2 ,C,H2 )
** Define fuel & oxidiser composition & temperatures
SCRS(FUIN,YFU1IN,YOX1IN,YH21I,YCO1I,YH2O1I,YCO21I,YN21IN,TEMFU)
SCRS(OXIN,YFU2IN,YOX2IN,YH22I,YCO2I,YH2O2I,YCO22I,YN22IN,TEMOX)
SCRS(PROP,CHEMKIN,SCRS)
ELSE
+ CHEMKIN(SPECIES,CH4,O2,H2,CO,H2O,CO2,H,O,OH,HO2,NO,N,N2)
+ NAME(45)=SRAD;NAME(44)=F
+ STORE(H1,SRAD,F)
ENDIF
MESG( Enter required combustion model: default mixed-is-burnt
MESG( The options are:
MESG( BURN - Mixed-is-burnt infinite-rate model
MESG( DDEL - Infinite-rate model with Double-Delta PDF
READVDU(CPDF,CHAR,BURN)
CASE :CPDF: OF
WHEN BURN,4
+ MESG(Mixed-is-burnt infinite rate model
WHEN DDEL,4
+ MESG(Infinite-rate model with Double-Delta PDF
+ SCRS(PDF,DDELTA)
ENDCASE
*** END OF EXTENDED SCRS MODEL SETTINGS
IF(BLOCK) THEN
+ CONPOR(BLK1,0.0,CELL,1,NX,NYG+1,NY,1,NZ)
+ OUTPUT(VPOR,P,P,Y,P,P,P)
ENDIF
TURMOD(KEMODL)
GROUP 8. Terms (in differential equations) & devices
GROUP 9. Properties of the medium (or media)
ENUL=3.E-5
IF(.NOT.RUN1) THEN
+ RHO1=GRND9;TMP1=0.0;TMP1B=13
+ STORE(KE,EP)
IF(BFC) THEN
+ NAME(42)=F;NAME(45)=SPH1;NAME(44)=YSUM;NAME(43)=SRAD
+ STORE(F,SRAD,YSUM,SPH1)
ENDIF
ENDIF
GROUP 11. Initialization of variable or porosity fields
INIADD=F
IF(RUN1) THEN
+ FIINIT(W1)=10.0
+ PATCH(INIT,INIVAL,1,NX,1,NYF,1,NZ,1,LSTEP)
+ INIT(INIT,W1,0.0,WINF); FIINIT(KE)=KEINF; FIINIT(EP)=EPINF
IF(THRAD) THEN
+ REAL(TGUESS);TGUESS=1791.0; FIINIT(SRAD)=0.07*SIGMA*TGUESS**4
ENDIF
IF(HSOLV) THEN
** inlet enthalpy as computed in GROUND
+ FIINIT(H1)=1.267E5
ENDIF
ELSE
+ RESTRT(ALL)
+ DO KK=CHSOB+6,CHSOB+11
+ FIINIT(:KK:)=1.E-6
+ ENDDO
ENDIF
GROUP 19. Data communicated by satellite to GROUND
(Moved to here so that it can be used in Group 13)
IF(.NOT.RUN1) THEN
+ CHEMKIN(REACT,TWOPNT)
+ CHSOC=1.0
+ CHEMKIN(PROP,RADN,CHSOC)
ENDIF
GROUP 13. Boundary conditions and special sources
** Fuel Stream Inlet Conditions
REAL(YIN);INTEGER(KK1)
IF(RUN1) THEN
INLET(SCRSF,LOW,1,NX,1,NYF,1,1,1,LSTEP)
+ VALUE(SCRSF,F,1.0)
+ VALUE(SCRSF,KE,KEINF); VALUE(SCRSF,EP,EPINF)
IF(HSOLV) THEN
+ VALUE(SCRSF,H1,GRND1)
ENDIF
IF(BFC) THEN
+ VALUE(SCRSF,P1,GRND3)
+ VALUE(SCRSF,V1,GRND3); VALUE(SCRSF,W1,GRND3)
+ VALUE(SCRSF,VCRT,ZERO); VALUE(SCRSF,WCRT,WINF)
ELSE
+ VALUE(SCRSF,P1,GRND1); VALUE(SCRSF,W1,WINF)
ENDIF
ELSE
+ INLET(NOCPCKF,LOW,1,NX,1,NYF,1,1,1,LSTEP)
+ VALUE(NOCPCKF,P1,GRND9); VALUE(NOCPCKF,W1,WINF)
+ VALUE(NOCPCKF,TMP1,TEMFU)
** inlet mass fractions
+ VALUE(NOCPCKF,NO,1.E-5)
+ VALUE(NOCPCKF,H2,7.283E-3); VALUE(NOCPCKF,CO,0.1035)
+ VALUE(NOCPCKF,H2O,0.1199); VALUE(NOCPCKF,CO2,0.06431)
ENDIF
** Oxidiser Stream Inlet Conditions
IF(RUN1) THEN
+ INLET(SCRSO,LOW,1,NX,NYF+1,NYG,1,1,1,LSTEP)
+ VALUE(SCRSO,P1,GRND1); VALUE(SCRSO,F,0.0)
+ VALUE(SCRSO,KE,KEINO); VALUE(SCRSO,EP,EPINO)
IF(HSOLV) THEN
+ VALUE(SCRSO,H1,GRND1)
ENDIF
IF(BFC) THEN
+ VALUE(SCRSO,P1,GRND3)
+ VALUE(SCRSO,V1,GRND3); VALUE(SCRSO,W1,GRND3)
+ VALUE(SCRSO,VCRT,ZERO); VALUE(SCRSO,WCRT,WINO)
ELSE
+ VALUE(SCRSO,P1,GRND1); VALUE(SCRSO,W1,WINO)
ENDIF
ELSE
+ INLET(NOCPCKO,LOW,1,NX,NYF+1,NYG,1,1,1,LSTEP)
+ VALUE(NOCPCKO,P1,GRND9); VALUE(NOCPCKO,W1,WINO)
+ VALUE(NOCPCKO,TMP1,TEMOX)
** inlet mass fractions
+ VALUE(NOCPCKO,O2,9.784E-02); VALUE(NOCPCKO,H2O,7.110E-02)
+ VALUE(NOCPCKO,CO2,8.606E-02); VALUE(NOCPCKO,N2,7.450E-01)
ENDIF
** note that wall type patch names must have WALL for characters
2 to 5 inclusive for ground-coded H1 convective b.c. to work
** Outer Wall Boundary Conditions
PATCH(NWALL1,NWALL,1,NX,NYG,NYG,1,NZ,1,LSTEP)
COVAL(NWALL1,W1,GRND2,0.0);COVAL(NWALL1,KE,GRND2,GRND2)
COVAL(NWALL1,EP,GRND2,GRND2)
** convective heat transfer wall function
COVAL(NWALL1,H1,GRND2,GRND); RG(50)=TWAL
IF(THRAD) THEN
+ PATCH(RWALLN,NORTH,1,NX,NYG,NYG,1,NZ,1,LSTEP)
+ COVAL(RWALLN,SRAD,EMISW/(2.0-EMISW),EMPW)
+ PATCH(RADINF,LOW,1,NX,1,NYF,1,1,1,LSTEP)
+ COVAL(RADINF,SRAD,0.5,EMPG*TEMFU**4)
+ PATCH(RADINO,LOW,1,NX,NYF+1,NYG,1,1,1,LSTEP)
+ COVAL(RADINO,SRAD,0.5,EMPG*TEMOX**4)
ENDIF
** Exit Boundary Conditions
OUTLET(OUT,HIGH,1,NX,1,NYG,NZ,NZ,1,LSTEP)
REAL(CMOUT);CMOUT=10.0
IF(RUN1) THEN
+ COVAL(OUT,P1,CMOUT,0.0)
+ VALUE(OUT,V1,0.0); VALUE(OUT,W1,0.0); VALUE(OUT,F,0.0)
ELSE
** convert CMOUT from s/m**3 to g*s/(kg*cm**3) for nox run
+ COVAL(OUT,P1,CMOUT*1.E-3,0.0)
+ DO KK=CHSOB,CHSOB+12
+ VALUE(OUT,:KK:,0.0)
+ ENDDO
ENDIF
KELIN=1.0
GROUP 15. Termination of sweeps
IF(.NOT.RUN1) THEN
+ DO KK=CHSOB,CHSOB+12
+ RESREF(:KK:)=1.E-12*RESREF(:KK:)
+ ENDIT(:KK:)=1.E-20
+ ENDDO
ENDIF
IF(RUN1) THEN
+ LSWEEP=200
ELSE
+ LSWEEP=80
ENDIF
GROUP 16. Termination of iterations
GROUP 17. Under-relaxation devices
REAL(RLXFAC); RLXFAC=2.0*ZWLAST/WINF/NZ
IF(RUN1) THEN
+ RELAX(V1,FALSDT,RLXFAC); RELAX(W1,FALSDT,RLXFAC)
+ RELAX(KE,LINRLX,0.4); RELAX(EP,LINRLX,0.4)
+ RELAX(DEN1,LINRLX,0.5); RELAX(F,LINRLX,0.8)
IF(HSOLV) THEN
+ RELAX(H1,FALSDT,1.0)
ENDIF
IF(THRAD) THEN
+ RELAX(SRAD,FALSDT,1.0)
ENDIF
ELSE
** LG(94)=T so as to convert density from kg/m**3 to g/cm**3 (by
multiplying by 1.E-3) for use in CHEMKIN on 1st sweep of NOX
run. However, on NOX restart runs you must set LG(94)=F
+ LG(94)=T; RELAX(DEN1,LINRLX,1.E-10)
+ DO KK=CHSOB,CHSOB+12
+ RELAX(:KK:,FALSDT,1.E-2)
+ ENDDO
ENDIF
GROUP 18. Limits on variables or increments to them
Prevent mass-fractions from straying outside the physical
limits (0,1), and prevent the temperature from falling below
0 C (chosen to be the reference temperature)
IF(RUN1) THEN
+ DO KK=CHSOB,CHSOB+6
+ VARMIN(:KK:)=0.0; VARMAX(:KK:)=1.0
+ ENDDO
ELSE
+ DO KK=CHSOB,CHSOB+12
+ VARMIN(:KK:)=-1.E-7; VARMAX(:KK:)=1.0+1.E-7
+ ENDDO
+ OUTPUT(SRAD,P,P,P,P,N,N); OUTPUT(SPH1,P,P,P,P,N,N)
+ OUTPUT(VIST,P,P,P,P,N,N); OUTPUT(F,P,P,P,P,N,N)
ENDIF
OUTPUT(CH4,N,N,P,N,N,N)
VARMIN(TMP1)=1.E-10; VARMIN(DEN1)=1.E-10
OUTPUT(TMP1,P,P,P,P,Y,Y); OUTPUT(DEN1,P,P,P,P,Y,Y)
GROUP 20. Preliminary print-out
ECHO=T
GROUP 21. Print-out of variables
NYPRIN=1;NZPRIN=1
GROUP 22. Spot-value print-out
TSTSWP=-1;IYMON=NYF+2;IZMON=NZ-1;ITABL=3;NPLT=10
GROUP 23. Field print-out and plot control
GROUP 24. Dumps for restarts
DISTIL=T
RUN1
IF(:CPDF:.EQ.BURN) THEN
+ EX(P1 )=2.800E+01; EX(V1 )=6.527E-01; EX(W1 )=3.330E+01
+ EX(KE )=9.627E+01; EX(EP )=1.521E+04; EX(H1 )=3.221E+05
+ EX(CH4 )=0.000E+00; EX(O2 )=2.489E-02; EX(H2 )=1.093E-03
+ EX(CO )=1.553E-02; EX(H2O )=1.093E-01; EX(CO2 )=1.210E-01
+ EX(N2 )=7.281E-01; EX(F )=4.214E-01; EX(SRAD)=8.462E+04
+ EX(YSUM)=1.000E+00; EX(SPH1)=1.482E+03; EX(TMP1)=2.024E+03
+ EX(DEN1)=1.637E-01; EX(VIST)=7.197E-02
ENDIF
IF(:CPDF:.EQ.DDEL) THEN
+ EX(P1 )=2.587E+01; EX(V1 )=6.697E-01; EX(W1 )=3.367E+01
+ EX(KE )=9.665E+01; EX(EP )=1.516E+04; EX(H1 )=1.958E+05
+ EX(CH4 )=0.000E+00; EX(O2 )=2.675E-02; EX(H2 )=1.202E-03
+ EX(CO )=1.708E-02; EX(H2O )=1.083E-01; EX(CO2 )=1.185E-01
+ EX(N2 )=7.282E-01; EX(F )=4.209E-01; EX(SRAD)=3.236E+04
+ EX(YSUM)=1.000E+00; EX(SPH1)=1.488E+03; EX(TMP1)=2.094E+03
+ EX(DEN1)=1.596E-01; EX(VIST)=7.255E-02
ENDIF
RUN2
IF(:CPDF:.EQ.BURN) THEN
+ EX(P1 )=2.621E+01; EX(V1 )=6.660E-01; EX(W1 )=3.400E+01
+ EX(KE )=9.693E+01; EX(EP )=1.530E+04; EX(H1 )=2.029E+05
+ EX(CH4 )=1.066E-20; EX(O2 )=2.681E-02; EX(H2 )=1.121E-03
+ EX(CO )=2.161E-02; EX(H2O )=1.077E-01; EX(CO2 )=1.118E-01
+ EX(H )=5.525E-05; EX(O )=2.730E-04; EX(OH )=2.164E-03
+ EX(HO2 )=1.667E-06; EX(NO )=8.744E-04; EX(N )=8.965E-09
+ EX(N2 )=7.277E-01; EX(SRAD)=3.282E+04; EX(SPH1)=1.490E+03
+ EX(TMP1)=2.110E+03; EX(DEN1)=1.588E-04; EX(YSUM)=1.000E+00
+ EX(VIST)=7.240E-02; EX(F)=4.235E-01
ENDIF
IF(RUN1) THEN
+ SELREF=T;NSAVE=COH2
ELSE
+ SELREF=F;NAMFI=COH2;NSAVE=RNOX; FIINIT(N2)=0.71
ENDIF
STORE(DFDY,DFDZ,GENF)
OUTPUT(YSUM,P,P,P,P,Y,Y)