C406: Combustion and NOX formation in a burner


 TALK=T;RUN( 1, 1)
  PHOTON USE
  p

  1 15 1


  vec z 1 sh
  gr ou z 1
  gr ou z 1 1 x 7 7 y 11 15
  gr ou z 1 1 x 7 7 y 4 8
  gr ou z 1 1 x 1 4 y 4 4
  msg( Velocity vectors
  pause;cl

  con h1 z 1 fil;.001
  gr ou z 1
  msg( Enthalpy contours
  pause;cl
  con tmp1 z 1 fil;.001
  gr ou z 1
  msg( Temperature
  pause;cl
  con den1 z 1 fil;.001
  gr ou z 1
  msg( Density
  pause;cl
  con ych4 z 1 fil;.001
  gr ou z 1
  msg( CH4 mass fraction
  pause;cl
  con yo2 z 1 fil;.001
  gr ou z 1
  msg( O2 mass fraction
  pause;cl
  con yco2 z 1 fil;.001
  gr ou z 1
  msg( CO2 mass fraction
  pause;cl
  con yco z 1 fil;.001
  gr ou z 1
  msg( CO mass fraction
  pause;cl
  con yh2o z 1 fil;.001
  gr ou z 1
  msg( H2O mass fraction
  pause;cl
  con yn2 z 1 fil;.001
  gr ou z 1
  msg(  N2 mass fraction
  pause;cl
  con ynox z 1 fil;.001
  gr ou z 1
  msg(  NOX mass fraction
  pause;cl
  con rnox z 1 fil;.001
  gr ou z 1
  msg( Volumetric rate of NOX formation
  pause;cl
  con rch z 1 fil;.001
  gr ou z 1
  msg( Volumetric rate of CH4 combustion
  pause;cl
  con rco z 1 fil;.001
  gr ou z 1
  msg( Volumetric rate of CO combustion
  pause;cl
  ENDUSE

  DISPLAY

  
  COMBUSTION and NOX FORMATION IN A BURNER

  *   Two-step reaction of combustion:

      2CH4  + 3O2 ->  2CO +4H2O
      2CO   +  O2 ->  2CO2

  *   Zeldovich's, partial-equilibrium mechanism of NOX formation

  ENDDIS
     GROUP 1. Run title and other preliminaries
     -------------------------------
TEXT(Two-step reaction: 1 fluid
     -------------------------------
REAL(HIN,GMIXL,CLEN,WIDTH,WIN1,WIN2,REYNO,WD2)
REAL(TKEIN1,EPIN1,TKEIN2,EPIN2)
INTEGER(IYJ);IYJ=3
REYNO=1.E6;WIDTH=0.3;HIN=1.;WD2=0.5*WIDTH
WIN1=4.;WIN2=4.

    GROUP 3. X-direction grid specification
NX=20;CLEN=70.*WD2;GRDPWR(X,NX,CLEN,2.0)
    GROUP 4. Y-direction grid specification
NY=15;GRDPWR(Y,NY,WD2,1.0)

    GROUP 5. Z-direction grid specification

    GROUP 7. Variables stored, solved & named

       * Solve for P1, U1, V1 and H1  - specific enthalpy
SOLVE(P1,U1,V1,H1)
SOLUTN(P1,Y,Y,Y,N,N,N)
       * Solve for single-fluid  mass fractions
SOLVE(YCH4,YO2,YCO,YH2O,YCO2,YN2,YNOX)
       * Store reciprocal of turbulent time scale
STORE(EPKE)
       * Combustion related inputs
REAL(CPFU,CPOX,HCH,HCO,HFUEL,TFUEL,TOX,HOX)
         ** Inlet temperatures of fuel and oxidant
TFUEL = 300.0;TOX   =600.0
         ** Varying specific heats
CPFU=  1059.+0.275*(TFUEL-300.)
CPOX=  1059.+0.275*(TOX-300.)
         ** Heat of combustion
HCH=5.5e7
HCO=7.45e6
         ** Inlet enthalpies of fuel and oxidant
HOX  = CPOX*TOX;HFUEL= CPFU*TFUEL
       * Select K-E model of hydrodynamic turbulencs
TURMOD(KEMODL)
KELIN=3.
    GROUP 8. Terms (in differential equations) & devices
TERMS(H1,N,Y,Y,N,Y,N)
    GROUP 9. Properties of the medium (or media)
REAL(RHOIN1,RHOIN2)
RHOIN1=1.e5*29./(8314*TFUEL)
RHOIN2=1.e5*29./(8314*TOX)
ENUL=WIN1*WIDTH/REYNO
store(den1)
 RHO1    = GRND5
 RHO1A   = 0.000000E+00 ;RHO1B  = 3.496503E-03
 RHO1C   = 7.142857E-01
 PRESS0  = 1.000000E+05

    Group 10.

    PLANTBEGIN
NAMSAT=MOSG
           ** Temperature
store(TMP1)
TMP1=GRND
   TEMP1=H1/(1059.+0.275*(TMP1-300.))
    * Reaction rate sources
           ** Reaction:   2CH4  + 3O2 ->  2CO +4H2O
PATCH(CH42CO,VOLUME,1,NX,1,NY,1,NZ,1,1)
  CO=Rch/(YCH4+tiny)
COVAL(CH42CO,YCH4,GRND,0.0)
  CO=3.*Rch/(YO2+tiny)
COVAL(CH42CO,YO2,GRND,0.0)
  VAL=1.75*Rch
COVAL(CH42CO,YCO,FIXFLU,GRND)
  VAL= 2.25*Rch
COVAL(CH42CO,YH2O,FIXFLU,GRND)
           ** Reaction:   2CO +  O2 -> 2CO2
PATCH(CO2CO2,VOLUME,1,NX,1,NY,1,NZ,1,1)
  CO=Rco/(YCO+tiny)
COVAL(CO2CO2,YCO,GRND,0.)
  CO=0.57*Rco/(YO2+tiny)
COVAL(CO2CO2,YO2,GRND,0.)
  VAL=1.57*Rco
COVAL(CO2CO2,YCO2,FIXFLU,GRND)
           ** Energy release with CO and CH4
PATCH(HEATCO2S,VOLUME,1,NX,1,NY,1,NZ,1,1)
   VAL=:Hco:*Rco+:Hch:*Rch
COVAL(HEATCO2S,H1,FIXFLU,GRND)
           ** NOX formation
PATCH(NOXSOR,VOLUME,1,NX,1,NY,1,NZ,1,1)
     VAL=RNOX
COVAL(NOXSOR,YNOX,FIXFLU,GRND)
           ** Combustion reaction rates
store(Rch,Rco,Rchk,Rcok)
   Rch=4.*EPKE*AMIN1(YCH4,YO2/3.)*DEN1
   Rco=4.*EPKE*AMIN1(YCO,YO2/.57)*DEN1
   Rchk=DEN1*1.15e9*EXP(-2.4e4/AMAX1(100.,TMP1))*$
               YO2**1.3/(YCH4**0.3+tiny)*DEN1
   Rcok=DEN1*5.4e9*EXP(-1.5e4/AMAX1(100.,TMP1))*$
              YCO*YO2**0.25*YH2O**0.5*DEN1
   Rch=amin1(Rch,Rchk)
   Rco=amin1(Rco,Rcok)
          ** NOX reaction rates
store(FK1,MOX,RNOX)
               *** Single-fluid reaction constant
   FK1=1.8e8*EXP(-38370./AMAX1(TMP1,100.))
               *** Mole O atom concentration
   MOX=3.97e5*(YO2*DEN1/32./AMAX1(TMP1,100.))**.5 $
              *EXP(-31090./AMAX1(TMP1,100.))
               *** NOX reaction rate
   RNOX=30.*2.*FK1*MOX*YN2/28.*DEN1
           ** Check gofr mass conservation
store(sums)
   SUMS=YCH4+YO2+YCO+YH2O+YCO2+YN2


    PLANTEND

    GROUP 11. Initialization of variable or porosity fields
FIINIT(U1)=0.5*(WIN1+WIN2)
FIINIT(H1)=1.e7

CONPOR(BLK1,0.0,CELL,-7,-7,4,8,1,1)
CONPOR(BLK2,0.0,CELL,-1,-4,-4,-4,1,1)
CONPOR(BLK3,0.0,CELL,-7,-7,11,15,1,1)

  ** TKEIN = 0.25*WIN1*WIN1*FRIC where FRIC=0.018 AT REYNO=1.E5
TKEIN1=0.25*WIN1*WIN1*0.018
TKEIN2=0.25*WIN2*WIN2*0.018
  ** EPIN = 0.1643*KIN**1.5/LMIX where LMIX=0.045*WIDTH
GMIXL=0.011*WD2
EPIN2=TKEIN2**1.5/GMIXL*0.1643
EPIN1=TKEIN1**1.5/GMIXL*0.1643
    GROUP 13. Boundary conditions and special sources
  ** Inlet Boundaries
INLET(IN1,WEST,1,1,1,IYJ,1,NZ,1,1)
VALUE(IN1,P1 , RHOIN1*WIN1)
VALUE(IN1,U1 , WIN1)
VALUE(IN1,H1,HFUEL)
VALUE(IN1,YCH4,1.0)
VALUE(IN1,YO2 ,0.0)
VALUE(IN1,YCO ,0.0)
VALUE(IN1,YH2O,0.0)
VALUE(IN1,YCO2,0.0)
VALUE(IN1,YN2 ,0.0)
VALUE(IN1,YNOX ,0.0)
VALUE(IN1,KE , TKEIN1)
VALUE(IN1,EP , EPIN1)

INLET(IN2,WEST,1,1,IYJ+1,NY,1,NZ,1,1)
VALUE(IN2,P1, RHOIN2*WIN2)
VALUE(IN2,U1, WIN2)
VALUE(IN2,H1,HOX)
VALUE(IN2,YCH4,0.0)
VALUE(IN2,YO2 ,0.232)
VALUE(IN2,YCO ,0.0)
VALUE(IN2,YH2O,0.0)
VALUE(IN2,YCO2,0.0)
VALUE(IN2,YN2 ,0.768)
VALUE(IN2,YNOX,0.0)
VALUE(IN2,KE, TKEIN2)
VALUE(IN2,EP, EPIN2)
   * Outlet boundary
PATCH(OUTLET,EAST,NX,NX,1,NY,1,NZ,1,1)
COVAL(OUTLET,P1,fixp,0.0)
   * North-Wall boundary (generalised wall functions)
WALL (WFNN,NORTH,1,NX,NY,NY,1,NZ,1,1)
    GROUP 15. Termination of sweeps
LSWEEP=500
RESFAC=1.e-3
    GROUP 16. Termination of iterations
LITHYD=10
 VARMIN(YCH4)=0.0;VARMAX(YCH4)=1.
 VARMIN(YO2) =0.0;VARMAX(YO2) =1.
 VARMIN(YCO) =0.0;VARMAX(YCO) =1.
 VARMIN(YH2O)=0.0;VARMAX(YH2O)=1.
 VARMIN(YCO2)=0.0;VARMAX(YCO2)=1.
 VARMIN(YN2) =0.0;VARMAX(YN2) =1.
 VARMIN(YNOX) =0.0;VARMAX(YNOX) =1.
 VARMIN(DEN1)=0.01;VARMAX(DEN1)=1.5
 VARMIN(H1)=1.e5;VARMAX(H1)=1.e8

    GROUP 17. Under-relaxation devices
RELAX(P1,LINRLX,.3)
RELAX(DEN1,LINRLX,0.03)
RELAX(V1,FALSDT,0.001)
RELAX(U1,FALSDT,0.001)
RELAX(KE,FALSDT,0.001)
RELAX(EP,FALSDT,0.001)
RELAX(H1,FALSDT,0.001)
RELAX(YCH4,FALSDT,0.001)
RELAX(YCO ,FALSDT,0.001)
RELAX(YCO2,FALSDT,0.001)
RELAX(YH2O,FALSDT,0.001)
RELAX(YCO2,FALSDT,0.001)
RELAX(YN2 ,FALSDT,0.001)
RELAX(YNOX ,FALSDT,0.001)

    GROUP 19. Data communicated by SATELLITE to GROUND

    GROUP 21. Print-out of variables
WALPRN=T;OUTPUT(KE,Y,Y,Y,Y,Y,Y)

    GROUP 22. Monitor print-out
IXMON=NX/2;IYMON=1;UWATCH=T

    GROUP 23. Field print-out and plot control
NPLT=1;NXPRIN=1;NYPRIN=1
NYPRIN=1;IYPRF=1;IYPRL=30
TSTSWP=-1
STOP