** LOAD(U111) from the USER Input Library
     ****** To load case: type LOAD(C111) ******
 
  PHOTON USE
  p
  PHI
  1 1 2
 
 
 
 
  use patgeo
  set con fi dep 3
  con mixf x 10 fi;0.001
  gr ou x 10
  te
  1
 
  mixture fraction contours at ix=10 plane
   0.50627E+03 0.31916E+03 CR
  msg Press  to continue
  pause
  text cl;con cl;redr
  con yco2 x 10 fi;0.001
  gr ou x 10
  te
  1
 
  CO2 contours at ix=10 plane
   0.50627E+03 0.31916E+03 CR
  msg Press  to continue
  pause
  text cl;con cl;redr
  con yco x 10 fi;0.001
  gr ou x 10
  te
  1
 
  CO contours at ix=10 plane
   0.50627E+03 0.31916E+03 CR
  msg Press  to continue
  pause
  text cl;con cl;redr
  con yo2 x 10 fi;0.001
  gr ou x 10
  te
  1
 
  O2 contours at ix=10 plane
   0.50627E+03 0.31916E+03 CR
  msg Press  to continue
  pause
  text cl;con cl;redr
  con tmp1 x 10 fi;0.001
  gr ou x 10
  te
  1
 
  gas phase temperature contours at ix=10 plane
   0.50627E+03 0.31916E+03 CR
  msg Press  to continue
  pause
  text cl;con cl;redr
  con tmp2 x 10 fi;0.001
  gr ou x 10
  te
  1
 
  solid phase temperature contours at ix=10 plane
   0.50627E+03 0.31916E+03 CR
  msg Press  to continue
  pause
  text cl;con cl;redr
  vec x 10 sh
  gr ou x 10
  te
  1
 
  velocity vectors at ix=10 plane
   0.50627E+03 0.31916E+03 CR
  msg Press  to continue
  pause
  enduse
    GROUP 1. Run title and other preliminaries
TEXT(3-D Coal Combustion In A Furnace  
TITLE
  DISPLAY
   The case considered is the 3d turbulent flow of combusting
   pulverised coal particles in a wall-fired power-station
   furnace. The calculation exploits two vertical planes of
   symmetry, so that only one quarter of the furnace is simulated.
   The burner arrangement comprises 3 rows of 4 horizontally-fired
   burners on opposing vertical walls, so that only 3 rows of 2 on
   a single wall are modelled in the actual solution domain.  Each
   burner consists of 3 air streams; tertiary air enters through a
   central air stream; swirling primary air carrying coal particles
   enters through an inner annulus located above the tertiary-air
   stream; and swirling secondary air enters through an outer
   annulus located above the primary air stream. The burner inflow
   boundary conditions are coded in the special ground routine
   furngr.f. The main combustion process is simulated using an
   IPSA-based 'equilibrium' 2-phase combustion model. Thermal
   radiation is accounted for by means of the radiosity model. A 
   2nd post-processing run can be performed which computes thermal 
   NOx based on CHO equilibrium as determined by the CREK module.
  ENDDIS
   The main combustion run requires 2000 sweeps for convergence
   and takes 2 hours on a 200MHz Pentium II. 
   ================================================
   Coal-combustion model; solid phase - carbon only
   reactions: C (s) + 0.5 O2  > CO      (exothermic)
              CO    + 0.5 O2  > CO2     (exothermic)
              C (s) + CO2     > 2CO     (endothermic)
              C (s) + H2O     > CO + H2 (endothermic)
              H2    + 0.5 O2  > H2O     (exothermic)
   ================================================
   noxdbegin
     ns   7
     nr   3
   noxdend
  ** control variables
BOOLEAN(OFA,NOXCAL,ANNINL,THRAD)
  ========================================
   Set appropriate switches
  ========================================
  ** Switch for overfire air ports
  -- set to "f" for no overfire air ports
  -- set to "t" for overfire air ports
OFA=F
  ** Switch for NOX calculation
  -- set to "f" for main combustion calculation,
  --        "t" for NOX post-processing calculation.
CHAR(NOXC)
MESG(
MESG(Main combustion calculation (M) or NOX Post-processing (P)
MESG(  (Default  =  M)
READVDU(NOXC,CHAR,M)
IF(:NOXC:.EQ.P)THEN
TEXT(NOX Post-processing calculation   : C111
+NOXCAL=T
ELSE
+NOXCAL=F
ENDIF
  ** Switch for annular inlet flow conditions
  -- set to "f" to use mass averaged per cell values of flow rates
  -- set to "t" to compute flow conditions in ground
ANNINL=T
   ** definition of integer variables
INTEGER(JXDBCK,JYDBCK,JZDBCK,JSDBCK,NS,NR,NBN)
INTEGER(IXB1,IXB2,IXB3,IXB4,IXB5,IXB6,IXFA1,IXFA2)
INTEGER(IYB1,IYB2,IYB3,IYB4,IYB5,IYB6)
INTEGER(IYFA1,IYFA2)
  ** definition of real variables
REAL (PI,AIN1,AIN2,AIN3,AOFA,CONV1,CONV2,CONV3,CONV4)
REAL (DTF,RIN1, RIN2, RIN3,ROFA, ROUT, ZQUARL)
REAL (FGIN1,FGIN2,FGIN3,FGOFA,FGOFAC,WGIN1,WGIN2,WGIN3,WGOFA)
REAL (FVIN1,FVIN2,FVIN3,FVOFA)
REAL (FCOAL,FMCOAL,WCOAL,ROCOAL,RMEAN,GABSR)
REAL (ALFC,ALFH2,ALFN2,TSIN,TGIN1,TGIN2,TGIN3,TGOFA)
REAL (HCOAL,HGIN1,HGIN2,HGIN3,HGOFA)
REAL (OMEGT,OMEGP,OMEGS,OMEGC,OMGOFA)
REAL (KEIN1,KEIN2,KEIN3,KEOFA,EPIN1,EPIN2,EPIN3,EPOFA)
REAL (DIST1,DIST2,DIST3,DISOFA)
REAL (BX1,BX2,HWID,HDEP,BY1,BYSPAC,YTOP)
REAL (FGINT,FGINC,FSINC,ASQ,WTOTIN,HTOTIN,KTOTIN,ETOTIN)
REAL (HGINAV,WGINAV,KGINAV,EGINAV)
REAL (XHWW,YHWW1,YHWW2,ZHWW1,ZHWW2,BNCN,TWALL,HWAL)
REAL (DTUBE,KWALL,PRWALL,MUWALL,OUTCOL)
REAL (XCB1,XCB2,XCB3,XCB4,XCB5,XCB6,XCOFA1,XCOFA2)
REAL (YCB1,YCB2,YCB3,YCB4,YCB5,YCB6,YCOFA1,YCOFA2)
       ======================================
  ***  Note that SI units are used throughout  ***
       ======================================
  ** furnace geometrical parameters  -  metres
  -- burner lateral positions
BX1=4.597; BX2=BX1+3.645
  -- height of lowest burners above furnace bottom
BY1=13.878
  -- burner vertical spacing
BYSPAC=3.048
  -- height of furnace top (as modelled) above top burners
YTOP=23.774
  -- half-width and half-depth of furnace
HWID=10.509; HDEP=6.398
  -- burners center coordinates
XCB1=BX1;YCB1=BY1
XCB2=BX1;YCB2=BY1+BYSPAC
XCB3=BX1;YCB3=BY1+2.*BYSPAC
XCB4=BX2;YCB4=BY1
XCB5=BX2;YCB5=BY1+BYSPAC
XCB6=BX2;YCB6=BY1+2.*BYSPAC
  -- overfire air center coordinates
XCOFA1=BX1;YCOFA1=BY1+3.*BYSPAC
XCOFA2=BX2;YCOFA2=BY1+3.*BYSPAC
  ** conversion factors
PI=3.14159
  -- inches to metres
CONV1=0.0254
  -- pounds to kilograms
CONV2=0.4536
  -- feet to metres
CONV3=CONV1*12.
  -- cu.ft to cu.m
CONV4=CONV3**3
  ** parameters for wall heat transfer calculation  -  SI Units
  -- tube outer diameter  -  inches > metres
DTUBE=2. * CONV1
  -- gas conductivity near walls
KWALL=3.365E-2
  -- Prandtl number near walls
PRWALL=0.688
  -- dynamic viscoity near walls
MUWALL=2.286E-5
  ** inlet radii  -  inches
  -- radius of innermost inlet tube
RIN1= 16.
  -- outer radius of middle inlet annulus
RIN2= 26.
  -- outer radius of outer secondary air annulus
RIN3= 53.
  -- overfire air port radius
ROFA= 53.
  -- quarl depth from inlet plane to furnace wall
ZQUARL= 11.313
  -- outer radius of 2-D solution domain  -  m
ROUT=1.
  -- radii, converted to SI units
RIN1= RIN1*CONV1/2.; RIN2=RIN2*CONV1/2.; RIN3=RIN3*CONV1/2.
ROFA= ROFA*CONV1/2.
ZQUARL=ZQUARL*CONV1
  ** geometry pertaining to  hanging water walls
  -- bottom of hww rel. to top of top burner/overfire air port
YHWW1=(3.419+10)*CONV3-RIN3-BYSPAC
  -- bottom of hww slot rel. to bottom of hww
YHWW2=(17-3.419)*CONV3
  -- width of hww gap
ZHWW1=2.089*CONV3
  -- full width of partial hww
ZHWW2=16.641*CONV3-ZHWW1
  -- lateral location of hww
XHWW=(BX2+BX1)/2.
  -- half-distance between burners
BNCN=(BX2-BX1)/2.
  ** inlet flow areas
AIN1=PI*RIN1*RIN1; AIN2=PI*RIN2*RIN2-AIN1;
AIN3=PI*RIN3*RIN3-AIN2-AIN1
AOFA=PI*ROFA*ROFA
  -- area of "square" inlet, i.e. simple approximation
ASQ=4.*RIN3*RIN3
  ** length scales for inlet epsilon
DIST1=(RIN2-RIN1)/10.; DIST2=(RIN3-RIN2)/10.; DIST3=RIN1/10.
DISOFA=ROFA/10.
  ** inlet mass flow rates  -  lb/s
     =====================
  -- primary air
FGIN1= 10.48
  -- secondary air
FGIN2= 44.47
  -- tertiary air
FGIN3= 7.97
  -- total flow rate (excluding over-fire air)
FGINT=FGIN1+FGIN2+FGIN3
  -- overfire air
FGOFA=4.447
  ** inlet parameters for coal
     =========================
  -- coal mass flow rate  -  lb/s
FCOAL= 5.075
  -- velocity of coal particles  -  ft/s
WCOAL= 70.06
  -- coal density
ROCOAL= 1350.0
  -- mass fraction of carbon in the solid
ALFC=0.93
  -- mass fraction of hydrogen in the solid
ALFH2=1.0-ALFC ; ALFN2=0.0
  -- convert velocity to SI
WCOAL=WCOAL*CONV3
  -- convert to mass flow per unit area in SI units
     remember - primary is ain2, secondary ain3, tertiary ain1.
FGINT=FGINT*CONV2
FGIN1=FGIN1*CONV2/AIN2; FGIN2=FGIN2*CONV2/AIN3
FGIN3=FGIN3*CONV2/AIN1; FMCOAL=FCOAL*CONV2; FCOAL=FMCOAL/AIN2
FGOFA=FGOFA*CONV2
  == NOTE: fmcoal is kg/s, fcoal is kg/(m**2.s), fgofa is in kg/s
  ** number of cells in each of x,y directions across burner
NBN=4
  ** mass inflow rate per cell
FGINC= FGINT/(NBN*NBN);FSINC= FMCOAL/(NBN*NBN)
FGOFAC= FGOFA/(NBN*NBN)
  ** inlet volume flows  -  ft**3/s
     ==================
  -- primary air
FVIN1=160.5
  -- secondary air
FVIN2=1071.5
  -- tertiary air
FVIN3=192.1
  -- overfire air
FVOFA=107.15
  -- convert to m**3/s
FVIN1=FVIN1*CONV4; FVIN2=FVIN2*CONV4;FVIN3=FVIN3*CONV4
FVOFA=FVOFA*CONV4
  -- mean inlet density  -  kg/m**3
RMEAN=FGINT*ASQ/(FVIN1+FVIN2+FVIN3)
  -- calculate inlet velocities  -  m/s
WGIN1=FVIN1/AIN2; WGIN2=FVIN2/AIN3; WGIN3=FVIN3/AIN1
WGOFA=FVOFA/AOFA
  -- total inlet axial momentum flow rate  -  kg/s * m/s
  -- note that the areas cancel  (fvin=win*ain)
WTOTIN= FVIN1*FGIN1 + FVIN2*FGIN2 + FVIN3*FGIN3
WGINAV= WTOTIN/FGINT
IF(OFA) THEN
+ WTOTIN = WTOTIN + WGOFA*FGOFA
ENDIF
  ** angular velocities at inlet  -  radians/sec
  ** solid-body rotation assumed
  -- coal
OMEGC= 10.
  -- primary air
OMEGP= 10.
  -- secondary air
OMEGS= 10.
  -- tertiary air
OMEGT= 10.
  -- overfire air
OMGOFA= 0.
  ** inlet temperatures  -  deg F
  -- coal
TSIN = 150.
  -- primary air
TGIN1= 150.
  -- secondary air
TGIN2= 500.
  -- tertiary air
TGIN3= 500.
  -- overfire air
TGOFA = 500.
  ** wall temperature (assumed constant)  -  deg F
TWALL=800.
  -- convert to deg K
TGIN1=273.+(TGIN1-32.)*5./9.
TGIN2=273.+(TGIN2-32.)*5./9.
TGIN3=273.+(TGIN3-32.)*5./9.
TGOFA=273.+(TGOFA-32.)*5./9.
TSIN =273.+(TSIN -32.)*5./9.
TWALL=273.+(TWALL-32.)*5./9.
  -- inlet conditions for turbulence - assume 5% intensity
KEIN1=0.0025*WGIN1*WGIN1;KEIN2=0.0025*WGIN2*WGIN2
KEIN3=0.0025*WGIN3*WGIN3;KEOFA=0.0025*WGOFA*WGOFA
EPIN1=0.1643*KEIN1**1.5/DIST1
EPIN2=0.1643*KEIN2**1.5/DIST2
EPIN3=0.1643*KEIN3**1.5/DIST3
EPOFA=0.1643*KEOFA**1.5/DISOFA
KTOTIN= KEIN1*FGIN1*AIN2 + KEIN2*FGIN2*AIN3 + KEIN3*FGIN3*AIN1
ETOTIN= EPIN1*FGIN1*AIN2 + EPIN2*FGIN2*AIN3 + EPIN3*FGIN3*AIN1
KGINAV= KTOTIN/FGINT;EGINAV= ETOTIN/FGINT
IF(OFA) THEN
+ KTOTIN = KTOTIN + KEOFA*FGOFA
+ ETOTIN = ETOTIN + EPOFA*FGOFA
ENDIF
  ** number of species in equilibrium dissociation
NS=7
  ** number of reactions in nox reactions
NR=3
  ===================================
  START OF GROUP-BY-GROUP Q1 SETTINGS
  ===================================
    GROUP 2. Transience; time-step specification
    GROUP 3. X-direction grid specification
IXB1=2;IXB2=2;IXB3=2;IXFA1=2
IXB4=5;IXB5=5;IXB6=5;IXFA2=5
NREGX=6
IREGX=1; GRDPWR(X, 4,BX1-RIN3, -1.2)
IREGX=2; GRDPWR(X, NBN, 2*RIN3,  1.)
IREGX=3; GRDPWR(X, 2, BNCN-RIN3, 1.)
IREGX=4; GRDPWR(X, 2, BNCN-RIN3, 1.)
IREGX=5; GRDPWR(X, NBN, 2*RIN3,  1.)
IREGX=6; GRDPWR(X, 2, HWID-(BX2+RIN3), 1.4)
    GROUP 4. Y-direction grid specification
IYB1=3;IYB2=5;IYB3=7;IYFA1=9
IYB4=3;IYB5=5;IYB6=7;IYFA2=9
NREGY=12
IREGY=1; GRDPWR(Y, 5, BY1-RIN3*2, -1.5)
IREGY=2; GRDPWR(Y, 1, RIN3, 1.)
IREGY=3; GRDPWR(Y, NBN, RIN3*2, 1.)
IREGY=4; GRDPWR(Y, 2, BYSPAC-RIN3*2, 1)
IREGY=5; GRDPWR(Y, NBN, RIN3*2, 1.)
IREGY=6; GRDPWR(Y, 2, BYSPAC-RIN3*2, 1)
IREGY=7; GRDPWR(Y, NBN, RIN3*2, 1.)
IREGY=8; GRDPWR(Y, 2, BYSPAC-RIN3*2, 1)
IREGY=9; GRDPWR(Y, NBN, RIN3*2, 1.)
IREGY=10; GRDPWR(Y, 2, YHWW1, 1.)
IREGY=11; GRDPWR(Y, 2, YHWW2, 1.)
IREGY=12; GRDPWR(Y, 5, YTOP-RIN3-YHWW1-YHWW2, 1.3)
    GROUP 5. Z-direction grid specification
NREGZ=4; IG(11)=0
IREGZ=1; GRDPWR(Z, 1, ZHWW1/2. ,1.)
IREGZ=2; GRDPWR(Z, 1, ZHWW1/2., 1.)
IREGZ=3; GRDPWR(Z, 4, ZHWW2, 1.1)
IREGZ=4; GRDPWR(Z, 2, HDEP-ZHWW1-ZHWW2, 1)
    GROUP 6. Body-fitted coordinates or grid distortion
    GROUP 7. Variables stored, solved & named
  ** Solve for one pressure, two velocities, the volume
     fractions of the two phases and the "shadow" volume
     fraction of the second (denser) phase.
ONEPHS=F
SOLVE(R1,R2,RS)
NAME(R1)=GAS;NAME(R2)=FUE;NAME(RS)=SHAD
SOLVE(U1,V1,W1,U2,V2,W2)
SOLUTN(P1,Y,Y,Y,P,P,P)
  ** Provide storage for inter-phase mass transfer.
STORE(MDOT,CFIP)
STORE(YO2,YCO,YCO2,YN2,YH2,YH2O)
  ** store temperature and density
STORE(TMP1,TMP2,RHO1,ENUT)
STORE(PRPS,VPOR,EPOR,NPOR,HPOR)
  ** enthalpy
SOLVE(H1,H2)
  ** Solve additionally for the miture fraction, i.e. the quantity
     of phase-2 material which has entered phase 1.
SOLVE(C1,C2);NAME(C1)=MIXF
  ** NOX solution
IF (NOXCAL) THEN
+ STORE(P1,U1,U2,V1,V2,W1,W2,GAS,FUE,SHAD,KE,EP,H1,H2,MIXF,C2)
+ SOLVE(C3,C5); STORE(C4,C6); NAME(C3)=XN; NAME(C5)=XNO
+ STORE(XO,XO2,XH,XOH,CRKT,NOSR,EQUI,DEGF)
+ CINT(XNO)=0.0; CINT(XN)=0.0; CINT(C4)=0.0; CINT(C6)=0.0
   ** Suppress phase-diffusion terms for NOX run; as these
      are undefined in XN & XNO equations if the volume 
      fractions are not solved.
+ TERMS(FUE,P,P,N,P,P,P);TERMS(GAS,P,P,N,P,P,P)
ENDIF
    GROUP 8. Terms (in differential equations) & devices
TERMS(H1,N,Y,Y,N,P,P);TERMS(H2,N,Y,Y,N,P,P)
IF(NOXCAL) THEN
+ THRAD=F
ELSE
+ MESG( Thermal radiation required ? (default=Y)
+ READVDU(ANS,CHAR,Y)
IF(:ANS:.EQ.Y) THEN
+ THRAD=T
ELSE
+ THRAD=F
ENDIF
IF(THRAD) THEN
+ REAL(ABSORB,SCAT,SIGMA,EMPW,EMISW,EMISG,EMPG)
+ ABSORB=0.1/XULAST;SCAT=ABSORB;EMISG=0.07
+ SIGMA=5.6697E-8;EMISW=1.0
+ EMPW=SIGMA*TWALL**4;EMPG=SIGMA*EMISG
+ RADIAT(RADI,ABSORB,SCAT,H1)
+ SOLUTN(SRAD,P,P,Y,P,P,P)
+ SOLUTN(H2,P,P,y,P,P,P)
ENDIF
ENDIF
    GROUP 9. Properties of the medium (or media)
  ** densities and temperatures
RHO1=7GASES; RHO2=ROCOAL; PRESS0=1.E5; TEMP0=0.
RHO1A=ALFC ; RHO1B=ALFH2
  -- take cpsolid=cpgas=1.e3
CP1=1.E3;CP2=CP1
REAL(ROGIN1,ROGIN2,ROGIN3)
       TMP2=LINH; TMP2A=0.; TMP2B=1./CP2
ROGIN1=PRESS0/(287.41*TGIN1)
ROGIN2=PRESS0/(287.41*TGIN2)
ROGIN3=PRESS0/(287.41*TGIN3)
  ** turbulence model
IF (.NOT.NOXCAL) THEN
+ TURMOD(KEMODL) ; KELIN=1
ELSE
+ STORE(KE,EP)
ENDIF
    GROUP 10. Inter-phase-transfer processes and properties
    ** Set constant interphase friction factor and activate
       the calculation of the interphase mass transfer by:
CFIPS=GRND1 ; CFIPC=1.E5
CFIPS=1.E6
CINT(H1)=1.0
CINT(H2)=CINT(H1)
  ** Note that grnd3 is a new mdot option, making the mass-
     transfer rate proportional to (cmdtc-mixf), where
     cmdtc stands for the saturation value of mixf, i.e. the
     largest value which can be attained as a result of mass
     transfer.  The coding is to be found in fn9899.
REAL(GFF,GFS)
GFS=0.232/(0.232+ALFC*16.0/12.0)
GFF=0.02577/(0.9146-0.5926*ALFC)
CMDOT=GRND3; CMDTA=1.E3; CMDTC=GFS
CINT(MIXF)=0.; CINT(C2)=0.
PHINT(MIXF)=1.0;PHINT(H1)=GRND7;PHINT(H2)=GRND7
    GROUP 11. Initialization of variable or porosity fields
FIINIT(GAS)=0.99999; FIINIT(FUE)=0.00001; FIINIT(SHAD)=0.00001
FIINIT(V1)=2.; FIINIT(V2)=FIINIT(V1); FIINIT(P1)=0.0
FIINIT(U1)=0.; FIINIT(U2)=0.;FIINIT(W1)=1.
FIINIT(W2)=FIINIT(W1);FIINIT(MDOT)=1.E-8
FIINIT(KE)=KEIN2; FIINIT(EP)=EPIN2
FIINIT(VPOR)=1.0; FIINIT(EPOR)=1.0
FIINIT(NPOR)=1.0; FIINIT(HPOR)=1.0
  ** settings for sloped hopper
REAL(GTHETA,POINTY,POINTZ);INTEGER(SLTYP1)
GTHETA=55.0;INIADD=F
POINTY=27.0307*CONV3
POINTZ=0.0; SLTYP1=4
PATCH(SLOPE1,INIVAL,1,NX,1,NY,#1,NZ,1,1)
COVAL(SLOPE1,VPOR,0.0,GRND);COVAL(SLOPE1,EPOR,0.0,GRND)
COVAL(SLOPE1,NPOR,0.0,GRND);COVAL(SLOPE1,HPOR,0.0,GRND)
COVAL(SLOPE1,PRPS,0.0,GRND)
 
REAL(HCCO2,HCCO,HHH2O,HCHX)
  ** HCCO2 = heat of combustion for  C  +     O2 -> CO2
     HCCO  =  "    "     "       "   C  + 0.5 O2 -> CO
     HHH2O =  "    "     "       "   H2 + 0.5 O2 -> H2O
        H = CP*T + HCHX*YCHX + HCOCO2*YCO * HH2*YH2
HCCO2=32.792E6; HCCO=9.208E6; HHH2O=120.9E6
HCHX=ALFC*HCCO2 + ALFH2*HHH2O
HCHX
 
HGIN1=CP1*(TGIN1-TEMP0);HGIN2=CP1*(TGIN2-TEMP0)
HGIN3=CP1*(TGIN3-TEMP0);HGOFA=CP1*(TGOFA-TEMP0)
 
HCOAL=CP2*(TSIN -TEMP0) + HCHX
HTOTIN= HGIN1*FGIN1*AIN2 + HGIN2*FGIN2*AIN3 + HGIN3*FGIN3*AIN1
HGINAV= HTOTIN/FGINT
 
HWAL=CP1*(TWALL-TEMP0)
 
FIINIT(H1)=HGIN2; FIINIT(H2)=HCOAL
IF(THRAD) THEN
+ FIINIT(SRAD)=EMPG
ENDIF
FIINIT(RHO1)=ROGIN2
  ** initially none of the second phase has entered the first
FIINIT(MIXF)=0.; FIINIT(C2)=1.
  ** restarts
IF (.NOT.NOXCAL) THEN
  == use the next line to activate or de-activate restarts
  == for the main combustion run.
    RESTRT(ALL)
ENDIF
IF (NOXCAL) THEN
+ RESTRT(P1,U1,U2,GAS,FUE,SHAD,H1,H2,MIXF,C2,TMP1,TMP2)
+ RESTRT(RHO1,RMIX,YH2O,YH2,YN2,YCO2,YO2,YCO,CFIP,MDOT)
+ RESTRT(U1,U2,V1,V2,W1,W2)
ENDIF
    GROUP 12. Unused
  GROUP 13. Boundary conditions and special sources
REAL(G71);G71=1.0; OUTCOL=1.E-4
   ** NOX patches
IF (NOXCAL) THEN
  ** NOX source
+ PATCH(NOXSR,FREEVL,1,NX,1,NY,1,NZ,1,LSTEP)
+ COVAL(NOXSR,XNO,FIXFLU,GRND)
+ COVAL(NOXSR,XN ,GRND ,GRND)
  ** NOX inlet patches
IF (.NOT.ANNINL) THEN
+ PATCH(NOXIN1,CELL,#IXB1,#IXB1,#IYB1,#IYB1,1,1,1,LSTEP)
+ COVAL(NOXIN1,XNO,FGINC ,ZERO)
+ COVAL(NOXIN1,XN ,FGINC ,ZERO)
+ PATCH(NOXIN2,CELL,#IXB2,#IXB2,#IYB2,#IYB2,1,1,1,LSTEP)
+ COVAL(NOXIN2,XNO,FGINC ,ZERO)
+ COVAL(NOXIN2,XN ,FGINC ,ZERO)
+ PATCH(NOXIN3,CELL,#IXB3,#IXB3,#IYB3,#IYB3,1,1,1,LSTEP)
+ COVAL(NOXIN3,XNO,FGINC ,ZERO)
+ COVAL(NOXIN3,XN ,FGINC ,ZERO)
+ PATCH(NOXIN4,CELL,#IXB4,#IXB4,#IYB4,#IYB4,1,1,1,LSTEP)
+ COVAL(NOXIN4,XNO,FGINC ,ZERO)
+ COVAL(NOXIN4,XN ,FGINC ,ZERO)
+ PATCH(NOXIN5,CELL,#IXB5,#IXB5,#IYB5,#IYB5,1,1,1,LSTEP)
+ COVAL(NOXIN5,XNO,FGINC ,ZERO)
+ COVAL(NOXIN5,XN ,FGINC ,ZERO)
+ PATCH(NOXIN6,CELL,#IXB6,#IXB6,#IYB6,#IYB6,1,1,1,LSTEP)
+ COVAL(NOXIN6,XNO,FGINC ,ZERO)
+ COVAL(NOXIN6,XN ,FGINC ,ZERO)
+ PATCH(NOXFA1,CELL,#IXFA1,#IXFA1,#IYFA1,#IYFA1,1,1,1,LSTEP)
+ COVAL(NOXFA1,XNO,FGOFAC ,ZERO)
+ COVAL(NOXFA1,XN ,FGOFAC ,ZERO)
+ PATCH(NOXFA2,CELL,#IXFA2,#IXFA2,#IYFA2,#IYFA2,1,1,1,LSTEP)
+ COVAL(NOXFA2,XNO,FGOFAC ,ZERO)
+ COVAL(NOXFA2,XN ,FGOFAC ,ZERO)
ELSE
    i.e. if (anninl) then
+ PATCH(NOXIN1,CELL,#IXB1,#IXB1,#IYB1,#IYB1,1,1,1,LSTEP)
+ COVAL(NOXIN1,XNO,GRND ,ZERO)
+ COVAL(NOXIN1,XN ,GRND ,ZERO)
+ PATCH(NOXIN2,CELL,#IXB2,#IXB2,#IYB2,#IYB2,1,1,1,LSTEP)
+ COVAL(NOXIN2,XNO,GRND ,ZERO)
+ COVAL(NOXIN2,XN ,GRND ,ZERO)
+ PATCH(NOXIN3,CELL,#IXB3,#IXB3,#IYB3,#IYB3,1,1,1,LSTEP)
+ COVAL(NOXIN3,XNO,GRND ,ZERO)
+ COVAL(NOXIN3,XN ,GRND ,ZERO)
+ PATCH(NOXIN4,CELL,#IXB4,#IXB4,#IYB4,#IYB4,1,1,1,LSTEP)
+ COVAL(NOXIN4,XNO,GRND ,ZERO)
+ COVAL(NOXIN4,XN ,GRND ,ZERO)
+ PATCH(NOXIN5,CELL,#IXB5,#IXB5,#IYB5,#IYB5,1,1,1,LSTEP)
+ COVAL(NOXIN5,XNO,GRND ,ZERO)
+ COVAL(NOXIN5,XN ,GRND ,ZERO)
+ PATCH(NOXIN6,CELL,#IXB6,#IXB6,#IYB6,#IYB6,1,1,1,LSTEP)
+ COVAL(NOXIN6,XNO,GRND ,ZERO)
+ COVAL(NOXIN6,XN ,GRND ,ZERO)
+ PATCH(NOXFA1,CELL,#IXFA1,#IXFA1,#IYFA1,#IYFA1,1,1,1,LSTEP)
+ COVAL(NOXFA1,XNO,GRND ,ZERO)
+ COVAL(NOXFA1,XN ,GRND ,ZERO)
+ PATCH(NOXFA2,CELL,#IXFA2,#IXFA2,#IYFA2,#IYFA2,1,1,1,LSTEP)
+ COVAL(NOXFA2,XNO,GRND ,ZERO)
+ COVAL(NOXFA2,XN ,GRND ,ZERO)
ENDIF
+ PATCH(OUTLET,NORTH,1,NX,NY,NY,#1,NZ,1,LSTEP)
+ COVAL(OUTLET, P1 ,OUTCOL,ZERO)
ELSE
    i.e. if (.not.noxcal) then
  ** sources for main calculation
IF (ANNINL) THEN
    ** burner 1
PATCH(BURN1,CELL,#IXB1,#IXB1,#IYB1,#IYB1,1,1,1,LSTEP)
COVAL(BURN1, P1,FIXFLU,GRND);COVAL(BURN1, P2,FIXFLU,GRND)
COVAL(BURN1, W1,ONLYMS,GRND);COVAL(BURN1, W2,ONLYMS,GRND)
COVAL(BURN1, U1,ONLYMS,GRND);COVAL(BURN1, U2,ONLYMS,GRND)
COVAL(BURN1, V1,ONLYMS,GRND);COVAL(BURN1, V2,ONLYMS,GRND)
COVAL(BURN1, H1,ONLYMS,GRND);COVAL(BURN1, H2,ONLYMS,GRND)
COVAL(BURN1, KE,ONLYMS,GRND);COVAL(BURN1, EP,ONLYMS,GRND)
COVAL(BURN1,MIXF,ONLYMS,ZERO);COVAL(BURN1, C2,ONLYMS,1.00)
    ** burner 2
PATCH(BURN2,CELL,#IXB2,#IXB2,#IYB2,#IYB2,1,1,1,LSTEP)
COVAL(BURN2, P1,FIXFLU,GRND);COVAL(BURN2, P2,FIXFLU,GRND)
COVAL(BURN2, W1,ONLYMS,GRND);COVAL(BURN2, W2,ONLYMS,GRND)
COVAL(BURN2, U1,ONLYMS,GRND);COVAL(BURN2, U2,ONLYMS,GRND)
COVAL(BURN2, V1,ONLYMS,GRND);COVAL(BURN2, V2,ONLYMS,GRND)
COVAL(BURN2, H1 ,ONLYMS,GRND);COVAL(BURN2, H2 ,ONLYMS,GRND)
COVAL(BURN2, KE ,ONLYMS,GRND);COVAL(BURN2, EP ,ONLYMS,GRND)
COVAL(BURN2,MIXF,ONLYMS,ZERO);COVAL(BURN2, C2 ,ONLYMS,1.00)
    ** burner 3
PATCH(BURN3,CELL,#IXB3,#IXB3,#IYB3,#IYB3,1,1,1,LSTEP)
COVAL(BURN3, P1 ,FIXFLU,GRND);COVAL(BURN3, P2 ,FIXFLU,GRND)
COVAL(BURN3, W1 ,ONLYMS,GRND);COVAL(BURN3, W2 ,ONLYMS,GRND)
COVAL(BURN3, U1 ,ONLYMS,GRND);COVAL(BURN3, U2 ,ONLYMS,GRND)
COVAL(BURN3, V1 ,ONLYMS,GRND);COVAL(BURN3, V2 ,ONLYMS,GRND)
COVAL(BURN3, H1 ,ONLYMS,GRND);COVAL(BURN3, H2 ,ONLYMS,GRND)
COVAL(BURN3, KE ,ONLYMS,GRND);COVAL(BURN3, EP ,ONLYMS,GRND)
COVAL(BURN3,MIXF,ONLYMS,ZERO);COVAL(BURN3, C2 ,ONLYMS,1.00)
    ** burner 4
PATCH(BURN4,CELL,#IXB4,#IXB4,#IYB4,#IYB4,1,1,1,LSTEP)
COVAL(BURN4, P1 ,FIXFLU,GRND);COVAL(BURN4, P2 ,FIXFLU,GRND)
COVAL(BURN4, W1 ,ONLYMS,GRND);COVAL(BURN4, W2 ,ONLYMS,GRND)
COVAL(BURN4, U1 ,ONLYMS,GRND);COVAL(BURN4, U2 ,ONLYMS,GRND)
COVAL(BURN4, V1 ,ONLYMS,GRND);COVAL(BURN4, V2 ,ONLYMS,GRND)
COVAL(BURN4, H1 ,ONLYMS,GRND);COVAL(BURN4, H2 ,ONLYMS,GRND)
COVAL(BURN4, KE ,ONLYMS,GRND);COVAL(BURN4, EP ,ONLYMS,GRND)
COVAL(BURN4,MIXF,ONLYMS,ZERO);COVAL(BURN4, C2 ,ONLYMS,1.00)
    ** burner 5
PATCH(BURN5,CELL,#IXB5,#IXB5,#IYB5,#IYB5,1,1,1,LSTEP)
COVAL(BURN5, P1 ,FIXFLU,GRND);COVAL(BURN5, P2 ,FIXFLU,GRND)
COVAL(BURN5, W1 ,ONLYMS,GRND);COVAL(BURN5, W2 ,ONLYMS,GRND)
COVAL(BURN5, U1 ,ONLYMS,GRND);COVAL(BURN5, U2 ,ONLYMS,GRND)
COVAL(BURN5, V1 ,ONLYMS,GRND);COVAL(BURN5, V2 ,ONLYMS,GRND)
COVAL(BURN5, H1 ,ONLYMS,GRND);COVAL(BURN5, H2 ,ONLYMS,GRND)
COVAL(BURN5, KE ,ONLYMS,GRND);COVAL(BURN5, EP ,ONLYMS,GRND)
COVAL(BURN5,MIXF,ONLYMS,ZERO);COVAL(BURN5, C2 ,ONLYMS,1.00)
    ** burner 6
PATCH(BURN6,CELL,#IXB6,#IXB6,#IYB6,#IYB6,1,1,1,LSTEP)
COVAL(BURN6, P1 ,FIXFLU,GRND);COVAL(BURN6, P2 ,FIXFLU,GRND)
COVAL(BURN6, W1 ,ONLYMS,GRND);COVAL(BURN6, W2 ,ONLYMS,GRND)
COVAL(BURN6, U1 ,ONLYMS,GRND);COVAL(BURN6, U2 ,ONLYMS,GRND)
COVAL(BURN6, V1 ,ONLYMS,GRND);COVAL(BURN6, V2 ,ONLYMS,GRND)
COVAL(BURN6, H1 ,ONLYMS,GRND);COVAL(BURN6, H2 ,ONLYMS,GRND)
COVAL(BURN6, KE ,ONLYMS,GRND);COVAL(BURN6, EP ,ONLYMS,GRND)
COVAL(BURN6,MIXF,ONLYMS,ZERO);COVAL(BURN6, C2 ,ONLYMS,1.00)
   ** overfire air port 1
PATCH(OFAIR1,CELL,#IXFA1,#IXFA1,#IYFA1,#IYFA1,1,1,1,LSTEP)
COVAL(OFAIR1, P1 ,FIXFLU,GRND);COVAL(OFAIR1, U1 ,ONLYMS,GRND)
COVAL(OFAIR1, V1 ,ONLYMS,GRND);COVAL(OFAIR1, W1 ,ONLYMS,WGOFA)
COVAL(OFAIR1, H1 ,ONLYMS,HGOFA);COVAL(OFAIR1, KE ,ONLYMS,KEOFA)
COVAL(OFAIR1, EP ,ONLYMS,EPOFA);COVAL(OFAIR1,MIXF,ONLYMS,ZERO)
COVAL(OFAIR1, C2 ,ONLYMS,1.00)
   ** overfire air port 2
PATCH(OFAIR2,CELL,#IXFA2,#IXFA2,#IYFA2,#IYFA2,1,1,1,LSTEP)
COVAL(OFAIR2, P1 ,FIXFLU,GRND);COVAL(OFAIR2, U1 ,ONLYMS,GRND)
COVAL(OFAIR2, V1 ,ONLYMS,GRND);COVAL(OFAIR2, W1 ,ONLYMS,WGOFA)
COVAL(OFAIR2, H1 ,ONLYMS,HGOFA);COVAL(OFAIR2, KE ,ONLYMS,KEOFA)
COVAL(OFAIR2, EP ,ONLYMS,EPOFA);COVAL(OFAIR2,MIXF,ONLYMS,ZERO)
COVAL(OFAIR2, C2 ,ONLYMS,1.00)
ELSE
    i.e. if (.not.anninl) then
   ** burner 1
PATCH(BURN1,CELL,#IXB1,#IXB1,#IYB1,#IYB1,1,1,1,LSTEP)
COVAL(BURN1, P1 ,FIXFLU,FGINC);COVAL(BURN1, P2 ,FIXFLU,FSINC)
COVAL(BURN1, W1 ,ONLYMS,WGINAV);COVAL(BURN1, W2 ,ONLYMS,WCOAL)
COVAL(BURN1, H1 ,ONLYMS,HGINAV);COVAL(BURN1, H2 ,ONLYMS,HCOAL)
COVAL(BURN1, KE ,ONLYMS,KGINAV);COVAL(BURN1, EP ,ONLYMS,EGINAV)
COVAL(BURN1,MIXF,ONLYMS,ZERO);COVAL(BURN1, C2 ,ONLYMS,1.00)
    ** burner 2
PATCH(BURN2,CELL,#IXB2,#IXB2,#IYB2,#IYB2,1,1,1,LSTEP)
COVAL(BURN2, P1 ,FIXFLU,FGINC);COVAL(BURN2, P2 ,FIXFLU,FSINC)
COVAL(BURN2, W1 ,ONLYMS,WGINAV);COVAL(BURN2, W2 ,ONLYMS,WCOAL)
COVAL(BURN2, H1 ,ONLYMS,HGINAV);COVAL(BURN2, H2 ,ONLYMS,HCOAL)
COVAL(BURN2, KE ,ONLYMS,KGINAV);COVAL(BURN2, EP ,ONLYMS,EGINAV)
COVAL(BURN2,MIXF,ONLYMS,ZERO);COVAL(BURN2, C2 ,ONLYMS,1.00)
    ** burner 3
PATCH(BURN3,CELL,#IXB3,#IXB3,#IYB3,#IYB3,1,1,1,LSTEP)
COVAL(BURN3, P1 ,FIXFLU,FGINC);COVAL(BURN3, P2 ,FIXFLU,FSINC)
COVAL(BURN3, W1 ,ONLYMS,WGINAV);COVAL(BURN3, W2 ,ONLYMS,WCOAL)
COVAL(BURN3, H1 ,ONLYMS,HGINAV);COVAL(BURN3, H2 ,ONLYMS,HCOAL)
COVAL(BURN3, KE ,ONLYMS,KGINAV);COVAL(BURN3, EP ,ONLYMS,EGINAV)
COVAL(BURN3,MIXF,ONLYMS,ZERO);COVAL(BURN3, C2 ,ONLYMS,1.00)
    ** burner 4
PATCH(BURN4,CELL,#IXB4,#IXB4,#IYB4,#IYB4,1,1,1,LSTEP)
COVAL(BURN4, P1 ,FIXFLU,FGINC);COVAL(BURN4, P2 ,FIXFLU,FSINC)
COVAL(BURN4, W1 ,ONLYMS,WGINAV);COVAL(BURN4, W2 ,ONLYMS,WCOAL)
COVAL(BURN4, H1 ,ONLYMS,HGINAV);COVAL(BURN4, H2 ,ONLYMS,HCOAL)
COVAL(BURN4, KE ,ONLYMS,KGINAV);COVAL(BURN4, EP ,ONLYMS,EGINAV)
COVAL(BURN4,MIXF,ONLYMS,ZERO);COVAL(BURN4, C2 ,ONLYMS,1.00)
    ** burner 5
PATCH(BURN5,CELL,#IXB5,#IXB5,#IYB5,#IYB5,1,1,1,LSTEP)
COVAL(BURN5, P1 ,FIXFLU,FGINC);COVAL(BURN5, P2 ,FIXFLU,FSINC)
COVAL(BURN5, W1 ,ONLYMS,WGINAV);COVAL(BURN5, W2 ,ONLYMS,WCOAL)
COVAL(BURN5, H1 ,ONLYMS,HGINAV);COVAL(BURN5, H2 ,ONLYMS,HCOAL)
COVAL(BURN5, KE ,ONLYMS,KGINAV);COVAL(BURN5, EP ,ONLYMS,EGINAV)
COVAL(BURN5,MIXF,ONLYMS,ZERO);COVAL(BURN5, C2 ,ONLYMS,1.00)
    ** burner 6
PATCH(BURN6,CELL,#IXB6,#IXB6,#IYB6,#IYB6,1,1,1,LSTEP)
COVAL(BURN6, P1 ,FIXFLU,FGINC);COVAL(BURN6, P2 ,FIXFLU,FSINC)
COVAL(BURN6, W1 ,ONLYMS,WGINAV);COVAL(BURN6, W2 ,ONLYMS,WCOAL)
COVAL(BURN6, H1 ,ONLYMS,HGINAV);COVAL(BURN6, H2 ,ONLYMS,HCOAL)
COVAL(BURN6, KE ,ONLYMS,KGINAV);COVAL(BURN6, EP ,ONLYMS,EGINAV)
COVAL(BURN6,MIXF,ONLYMS,ZERO);COVAL(BURN6, C2 ,ONLYMS,1.00)
   ** overfire air port 1
PATCH(OFAIR1,CELL,#IXFA1,#IXFA1,#IYFA1,#IYFA1,1,1,1,LSTEP)
COVAL(OFAIR1, P1 ,FIXFLU,FGOFAC);COVAL(OFAIR1, W1 ,ONLYMS,WGOFA)
COVAL(OFAIR1, H1 ,ONLYMS,HGOFA);COVAL(OFAIR1, KE ,ONLYMS,KEOFA)
COVAL(OFAIR1, EP ,ONLYMS,EPOFA)
   ** overfire air port 2
PATCH(OFAIR2,CELL,#IXFA2,#IXFA2,#IYFA2,#IYFA2,1,1,1,LSTEP)
COVAL(OFAIR2, P1 ,FIXFLU,FGOFAC);COVAL(OFAIR2, W1 ,ONLYMS,WGOFA)
COVAL(OFAIR2, H1 ,ONLYMS,HGOFA);COVAL(OFAIR2, KE ,ONLYMS,KEOFA)
COVAL(OFAIR2, EP ,ONLYMS,EPOFA)
ENDIF
   ** if no overfire air ports set appropriate patches to skip
IF(.NOT.OFA) THEN
IF(NOXCAL) THEN
NOXFA1=SKIPA;NOXFA2=SKIPB
ENDIF
OFAIR1=SKIPC;OFAIR2=SKIPD
ENDIF
  ** wall roughness height - metres
WALLA=0.005
  ** west wall - friction
WALL (WALLW,WEST,1,1,1,NY,#1,NZ,1,LSTEP)
  ** low  wall - friction
WALL (WALLA, LOW,#01,#06,#01,#02,$1,$1,1,LSTEP)
WALL (WALLB, LOW,#01,#01,#03,#03,$1,$1,1,LSTEP)
WALL (WALLC, LOW,#03,#04,#03,#03,$1,$1,1,LSTEP)
WALL (WALLD, LOW,#06,#06,#03,#03,$1,$1,1,LSTEP)
WALL (WALLE, LOW,#01,#06,#04,#04,$1,$1,1,LSTEP)
WALL (WALLF, LOW,#01,#01,#05,#05,$1,$1,1,LSTEP)
WALL (WALLG, LOW,#03,#04,#05,#05,$1,$1,1,LSTEP)
WALL (WALLH, LOW,#06,#06,#05,#05,$1,$1,1,LSTEP)
WALL (WALLI, LOW,#01,#06,#06,#06,$1,$1,1,LSTEP)
WALL (WALLJ, LOW,#01,#01,#07,#07,$1,$1,1,LSTEP)
WALL (WALLK, LOW,#03,#04,#07,#07,$1,$1,1,LSTEP)
WALL (WALLL, LOW,#06,#06,#07,#07,$1,$1,1,LSTEP)
IF (.NOT.OFA) THEN
WALL (WALLM, LOW,#01,#06,#08,#12,$1,$1,1,LSTEP)
ELSE
WALL (WALLM, LOW,#01,#06,#08,#08,$1,$1,1,LSTEP)
WALL (WALLN, LOW,#01,#01,#09,#09,$1,$1,1,LSTEP)
WALL (WALLO, LOW,#03,#04,#09,#09,$1,$1,1,LSTEP)
WALL (WALLP, LOW,#06,#06,#09,#09,$1,$1,1,LSTEP)
WALL (WALLQ, LOW,#01,#06,#10,#12,$1,$1,1,LSTEP)
ENDIF
  ** low  wall - heat transfer
PATCH(HTLWLA, LOW,#01,#06,#01,#02,$1,$1,1,LSTEP)
PATCH(HTLWLB, LOW,#01,#01,#03,#03,$1,$1,1,LSTEP)
PATCH(HTLWLC, LOW,#03,#04,#03,#03,$1,$1,1,LSTEP)
PATCH(HTLWLD, LOW,#06,#06,#03,#03,$1,$1,1,LSTEP)
PATCH(HTLWLE, LOW,#01,#06,#04,#04,$1,$1,1,LSTEP)
PATCH(HTLWLF, LOW,#01,#01,#05,#05,$1,$1,1,LSTEP)
PATCH(HTLWLG, LOW,#03,#04,#05,#05,$1,$1,1,LSTEP)
PATCH(HTLWLH, LOW,#06,#06,#05,#05,$1,$1,1,LSTEP)
PATCH(HTLWLI, LOW,#01,#06,#06,#06,$1,$1,1,LSTEP)
PATCH(HTLWLJ, LOW,#01,#01,#07,#07,$1,$1,1,LSTEP)
PATCH(HTLWLK, LOW,#03,#04,#07,#07,$1,$1,1,LSTEP)
PATCH(HTLWLL, LOW,#06,#06,#07,#07,$1,$1,1,LSTEP)
IF (.NOT.OFA) THEN
PATCH(HTLWLM, LOW,#01,#06,#08,#12,$1,$1,1,LSTEP)
ELSE
PATCH(HTLWLM, LOW,#01,#06,#08,#08,$1,$1,1,LSTEP)
PATCH(HTLWLN, LOW,#01,#01,#09,#09,$1,$1,1,LSTEP)
PATCH(HTLWLO, LOW,#03,#04,#09,#09,$1,$1,1,LSTEP)
PATCH(HTLWLP, LOW,#06,#06,#09,#09,$1,$1,1,LSTEP)
PATCH(HTLWLQ, LOW,#01,#06,#10,#12,$1,$1,1,LSTEP)
ENDIF
 
IF(THRAD) THEN
+ REAL(CORA);CORA=2.*EMISW/(2.-EMISW)
+ COVAL(HTLWLA,SRAD,CORA,EMPW);COVAL(HTLWLA,SRAD,CORA,EMPW)
+ COVAL(HTLWLB,SRAD,CORA,EMPW);COVAL(HTLWLC,SRAD,CORA,EMPW)
+ COVAL(HTLWLD,SRAD,CORA,EMPW);COVAL(HTLWLE,SRAD,CORA,EMPW)
+ COVAL(HTLWLF,SRAD,CORA,EMPW);COVAL(HTLWLG,SRAD,CORA,EMPW)
+ COVAL(HTLWLH,SRAD,CORA,EMPW);COVAL(HTLWLI,SRAD,CORA,EMPW)
+ COVAL(HTLWLJ,SRAD,CORA,EMPW);COVAL(HTLWLK,SRAD,CORA,EMPW)
+ COVAL(HTLWLL,SRAD,CORA,EMPW);COVAL(HTLWLM,SRAD,CORA,EMPW)
ELSE
+ COVAL(HTLWLA,H1,GRND,HWAL);COVAL(HTLWLA,H1,GRND,HWAL)
+ COVAL(HTLWLB,H1,GRND,HWAL);COVAL(HTLWLC,H1,GRND,HWAL)
+ COVAL(HTLWLD,H1,GRND,HWAL);COVAL(HTLWLE,H1,GRND,HWAL)
+ COVAL(HTLWLF,H1,GRND,HWAL);COVAL(HTLWLG,H1,GRND,HWAL)
+ COVAL(HTLWLH,H1,GRND,HWAL);COVAL(HTLWLI,H1,GRND,HWAL)
+ COVAL(HTLWLJ,H1,GRND,HWAL);COVAL(HTLWLK,H1,GRND,HWAL)
+ COVAL(HTLWLL,H1,GRND,HWAL);COVAL(HTLWLM,H1,GRND,HWAL)
ENDIF
IF (OFA) THEN
IF(THRAD) THEN
+ COVAL(HTLWLN,SRAD,CORA,EMPW);COVAL(HTLWLO,SRAD,CORA,EMPW)
+ COVAL(HTLWLP,SRAD,CORA,EMPW);COVAL(HTLWLQ,SRAD,CORA,EMPW)
ELSE
+ COVAL(HTLWLN,H1,GRND,HWAL);COVAL(HTLWLO,H1,GRND,HWAL)
+ COVAL(HTLWLP,H1,GRND,HWAL);COVAL(HTLWLQ,H1,GRND,HWAL)
ENDIF
ENDIF
  ** hanging water walls, friction
  -- note this does both sides of partial wall
WALL(HWWA,EAST,%3,%3,#11,#11,#1,#3,1,LSTEP); HWWA=SKIPA
WALL(HWWB,EAST,%3,%3,#12,#12,#3,#3,1,LSTEP); HWWB=SKIPB
WALL(HWWC,EAST,NX,NX,#11,#11,#1,NZ,1,LSTEP)
WALL(HWWD,EAST,NX,NX,#12,#12,#3,NZ,1,LSTEP)
  ** hanging water wall, resistance to through flow
PATCH(HWWUA,EAST,%3,%3,#11,#11,#1,#3,1,LSTEP)
PATCH(HWWUB,EAST,%3,%3,#12,#12,#3,#3,1,LSTEP)
COVAL(HWWUA,U1,GRND,0.)
COVAL(HWWUB,U1,GRND,0.)
  ** heat transfer for hanging water walls, and west face
PATCH(HTEHWA,EAST,%3,%3,#11,#11,#1,#3,1,LSTEP)
PATCH(HTEHWB,EAST,%3,%3,#12,#12,#3,#3,1,LSTEP)
PATCH(HTEHWC,WEST,$4,$4,#11,#11,#1,#3,1,LSTEP)
PATCH(HTEHWD,WEST,$4,$4,#12,#12,#3,#3,1,LSTEP)
PATCH(HTEHWE,EAST,NX,NX,#11,#11,#1,NZ,1,LSTEP)
PATCH(HTEHWF,EAST,NX,NX,#12,#12,#3,NZ,1,LSTEP)
PATCH(HTEWST,WEST, 1, 1,1, NY, #1,NZ,1,LSTEP)
IF(THRAD) THEN
+ COVAL(HTEHWA,SRAD,CORA,EMPW);COVAL(HTEHWB,SRAD,CORA,EMPW)
+ COVAL(HTEHWC,SRAD,CORA,EMPW);COVAL(HTEHWD,SRAD,CORA,EMPW)
+ COVAL(HTEHWE,SRAD,CORA,EMPW);COVAL(HTEHWF,SRAD,CORA,EMPW)
+ COVAL(HTEWST,SRAD,CORA,EMPW)
ELSE
+ COVAL(HTEHWA,H1,GRND,HWAL);COVAL(HTEHWB,H1,GRND,HWAL)
+ COVAL(HTEHWC,H1,GRND,HWAL);COVAL(HTEHWD,H1,GRND,HWAL)
+ COVAL(HTEHWE,H1,GRND,HWAL);COVAL(HTEHWF,H1,GRND,HWAL)
+ COVAL(HTEWST,H1,GRND,HWAL)
ENDIF
  ** Outlet at north end
PATCH(OUTLET,NORTH,1,NX,NY,NY,#1,NZ,1,LSTEP)
COVAL(OUTLET, P1 ,OUTCOL,ZERO);COVAL(OUTLET, P2 ,OUTCOL*RHO2,ZERO)
COVAL(OUTLET, U1 ,ONLYMS,ZERO);COVAL(OUTLET, U2 ,ONLYMS,ZERO)
COVAL(OUTLET, V1 ,ONLYMS,ZERO);COVAL(OUTLET, V2 ,ONLYMS,ZERO)
COVAL(OUTLET, W1 ,ONLYMS,ZERO);COVAL(OUTLET, W2 ,ONLYMS,ZERO)
COVAL(OUTLET, H1 ,ONLYMS,SAME);COVAL(OUTLET, H2 ,ONLYMS,SAME)
COVAL(OUTLET, KE ,ONLYMS,SAME);COVAL(OUTLET, EP ,ONLYMS,SAME)
COVAL(OUTLET,MIXF,ONLYMS,SAME);COVAL(OUTLET, C2 ,ONLYMS,1.00)
PATCH(GRAVTY,PHASEM,1,NX,1,NY,1,NZ,1,LSTEP)
COVAL(GRAVTY,V1,FIXFLU,-9.81);COVAL(GRAVTY,V2,FIXFLU,-9.81)
ENDIF
    (end of if(.not.noxcal) )
    GROUP 14. Downstream pressure for PARAB=.TRUE.
    GROUP 15. Termination of sweeps
LSWEEP=2000
IF (NOXCAL) THEN
+ LSWEEP=400
ENDIF
LITER(P1)=30
LITER(U1)=4;LITER(U2)=4;LITER(V1)=4
LITER(V2)=2;LITER(W1)=3;LITER(W2)=3
LITER(FUE)=2;LITER(GAS)=2;LITER(SHAD)=2
    GROUP 16. Termination of iterations
  -- recover inlet mass flow rates in kg/s (not per unit area)
FGIN1=FGIN1*AIN2; FGIN2=FGIN2*AIN3; FGIN3=FGIN3*AIN1
  GROUP 17. Under-relaxation device
  ** relaxations on volume fraction
RELAX(SHAD,LINRLX,0.3)
RELAX(FUE,LINRLX,0.3);RELAX(GAS,LINRLX,0.3)
  ** relaxations on pressure, density etc
RELAX(RHO1,LINRLX,0.3);RELAX(MDOT,LINRLX,0.7)
  ** relaxations on other variables
DTF= HDEP/(NZ*WGIN1)
RELAX( W1 ,FALSDT,DTF);RELAX( W2 ,FALSDT,DTF)
RELAX( V1 ,FALSDT,DTF);RELAX( V2 ,FALSDT,DTF)
RELAX( U1 ,FALSDT,DTF);RELAX( U2 ,FALSDT,DTF)
RELAX( H1 ,FALSDT,2.*DTF);RELAX( H2 ,FALSDT,2.*DTF)
RELAX( KE ,FALSDT,DTF/5.);RELAX( EP ,FALSDT,DTF/5.)
RELAX(MIXF,FALSDT,DTF/2.)
IF(THRAD) THEN
+ RELAX(SRAD,FALSDT,10.*DTF)
ENDIF
    GROUP 18. Limits on variables or increments to them
VARMAX(MIXF)=CMDTC; VARMIN(MIXF)=0.
VARMAX( C2 )=1.   ; VARMIN( C2 )=0.
VARMAX(SHAD)=1.   ; VARMIN(SHAD)=1.E-10
VARMAX( FUE)=1.   ; VARMIN( FUE)=1.E-10
VARMAX( GAS)=1.   ; VARMIN( GAS)=1.E-10
VARMIN(RHO1)=1.E-10;VARMIN(H2)=HCOAL
    ** additional settings for NOX calculation
IF (NOXCAL) THEN
+ RELAX(XNO,FALSDT,5.E-2)
+ RELAX(XN ,FALSDT,5.E-2)
+ VARMIN(XNO)=0.0;VARMIN(XN)=0.0
+ VARMAX(XNO)=1.0;VARMAX(XN)=1.0
  -- low reference for NOX variables to suppress early cutoff
+ SELREF=F
+ RESREF(XNO)=RESREF(GAS) * 1.E-10
+ RESREF(XN )=RESREF(XNO)
ENDIF
    GROUP 19. Data communicated by satellite to GROUND
  **  LSG23 removed as LSG62 now does what is needed.
CSG9=COAL
LSG62=NOXCAL
NAMGRD=FURN
  ** values to be communicated to Ground
SPEDAT(SET,FURN,XCB1,R,XCB1);SPEDAT(SET,FURN,YCB1,R,YCB1)
SPEDAT(SET,FURN,XCB2,R,XCB2);SPEDAT(SET,FURN,YCB2,R,YCB2)
SPEDAT(SET,FURN,XCB3,R,XCB3);SPEDAT(SET,FURN,YCB3,R,YCB3)
SPEDAT(SET,FURN,XCB4,R,XCB4);SPEDAT(SET,FURN,YCB4,R,YCB4)
SPEDAT(SET,FURN,XCB5,R,XCB5);SPEDAT(SET,FURN,YCB5,R,YCB5)
SPEDAT(SET,FURN,XCB6,R,XCB6);SPEDAT(SET,FURN,YCB6,R,YCB6)
SPEDAT(SET,FURN,FGIN1,R,FGIN3/AIN1)
SPEDAT(SET,FURN,FGIN2,R,FGIN1/AIN2)
SPEDAT(SET,FURN,FGIN3,R,FGIN2/AIN3)
SPEDAT(SET,FURN,FCOAL,R,FCOAL);SPEDAT(SET,FURN,WGIN1,R,WGIN3)
SPEDAT(SET,FURN,WGIN2,R,WGIN1);SPEDAT(SET,FURN,WGIN3,R,WGIN2)
SPEDAT(SET,FURN,WCOAL,R,WCOAL);SPEDAT(SET,FURN,HGIN1,R,HGIN3)
SPEDAT(SET,FURN,HGIN2,R,HGIN1);SPEDAT(SET,FURN,HGIN3,R,HGIN2)
SPEDAT(SET,FURN,HCOAL,R,HCOAL);SPEDAT(SET,FURN,KEIN1,R,KEIN3)
SPEDAT(SET,FURN,KEIN2,R,KEIN1);SPEDAT(SET,FURN,KEIN3,R,KEIN2)
SPEDAT(SET,FURN,EPIN1,R,EPIN3);SPEDAT(SET,FURN,EPIN2,R,EPIN1)
SPEDAT(SET,FURN,EPIN3,R,EPIN2)
SPEDAT(SET,FURN,BX1,R,BX1);SPEDAT(SET,FURN,BX2,R,BX2)
SPEDAT(SET,FURN,BY1,R,BY1);SPEDAT(SET,FURN,BY1A,R,BY1+BYSPAC)
SPEDAT(SET,FURN,BY1B,R,BY1+2*BYSPAC)
SPEDAT(SET,FURN,DIAT,R,2.0*RIN1);SPEDAT(SET,FURN,DIAP,R,2.0*RIN2)
SPEDAT(SET,FURN,DIAS,R,2.0*RIN3);SPEDAT(SET,FURN,GTHETA,R,GTHETA)
SPEDAT(SET,FURN,POINTY,R,POINTY);SPEDAT(SET,FURN,POINTZ,R,POINTZ)
SPEDAT(SET,FURN,OMEGT,R,OMEGT);SPEDAT(SET,FURN,OMEGP,R,OMEGP)
SPEDAT(SET,FURN,OMEGS,R,OMEGS);SPEDAT(SET,FURN,OMEGC,R,OMEGC)
SPEDAT(SET,FURN,TWALL,R,TWALL);SPEDAT(SET,FURN,OUTCOL,R,OUTCOL)
SPEDAT(SET,FURN,OMGOFA,R,OMGOFA)
SPEDAT(SET,FURN,XCOFA1,R,XCOFA1);SPEDAT(SET,FURN,YCOFA1,R,YCOFA1)
SPEDAT(SET,FURN,XCOFA2,R,XCOFA2);SPEDAT(SET,FURN,YCOFA2,R,YCOFA2)
SPEDAT(SET,FURN,GROFA,R,FGOFA/AOFA);SPEDAT(SET,FURN,GWOFA,R,WGOFA)
SPEDAT(SET,FURN,GKEOFA,R,KEOFA);SPEDAT(SET,FURN,GEPOFA,R,EPOFA)
SPEDAT(SET,FURN,GHOFA,R,HGOFA)
SPEDAT(SET,FURN,DTUBE,R,DTUBE);SPEDAT(SET,FURN,KWALL,R,KWALL)
SPEDAT(SET,FURN,PRWALL,R,PRWALL);SPEDAT(SET,FURN,MUWALL,R,MUWALL)
SPEDAT(SET,FURN,G71,R,G71)
SPEDAT(SET,FURN,JBCOPT,I,SLTYP1)
SPEDAT(SET,FURN,ANNINL,L,ANNINL)
    GROUP 20. Preliminary print-out
    GROUP 21. Print-out of variables
OUTPUT( P1 ,Y,Y,Y,Y,Y,Y);OUTPUT( C2 ,N,N,P,N,N,N)
OUTPUT(MDOT,N,N,P,N,N,N);OUTPUT(CFIP,N,N,P,N,N,N)
OUTPUT(ENUT,N,N,P,N,N,N);OUTPUT(RHO1,N,N,P,N,N,N)
IF (NOXCAL) THEN
+ OUTPUT( C4 ,N,N,P,N,N,N);OUTPUT( C6 ,N,N,P,N,N,N)
+ OUTPUT( H1 ,N,N,P,N,N,N);OUTPUT( H2 ,N,N,P,N,N,N)
+ OUTPUT( U2 ,N,N,P,N,N,N);OUTPUT(SHAD,N,N,P,N,N,N)
+ OUTPUT( FUE,N,N,P,N,N,N);OUTPUT(TMP2,N,N,P,N,N,N)
+ OUTPUT(CFIP,N,N,P,N,N,N);OUTPUT(TMP2,N,N,P,N,N,N)
+ OUTPUT(TMP1,P,P,P,Y,Y,P);OUTPUT(XO2,P,P,P,Y,Y,P)
ENDIF
    GROUP 22. Spot-value print-out
NPRMON=LSWEEP;IXMON=4;IYMON=11;IZMON=7
    ** controls for CREK debug printout
   DBGNOX=F
   JXDBCK=IXMON; JYDBCK=IYMON; JZDBCK=IZMON; JSDBCK=LSWEEP-1
    GROUP 23. Field print-out and plot control
TSTSWP=10; NPLT=50; ITABL=3; ORSIZ=0.4; IPLTF=2
IF (NOXCAL) THEN
+ NPLT=2
ENDIF
TSTSWP=-1
    GROUP 24. Dumps for restarts
DISTIL=T
EX(P1)=2.634E+04;EX(U1)=3.091E+00;EX(U2)=3.050E+00
EX(V1)=7.427E+00;EX(V2)=7.334E+00;EX(W1)=6.170E+00
EX(W2)=6.228E+00;EX(GAS)=9.797E-01;EX(KE)=1.300E+01
EX(EP)=4.428E+01;EX(H1)=1.406E+06;EX(H2)=3.954E+07
EX(MIXF)=2.949E-02;EX(HPOR)=9.706E-01
EX(NPOR)=9.766E-01;EX(EPOR)=9.663E-01;EX(VPOR)=9.663E-01
EX(TMP2)=1.401E+03;EX(TMP1)=1.406E+03;EX(YH2O)=1.858E-02
EX(YN2)=7.298E-01;EX(YCO2)=1.006E-01;EX(YO2)=1.308E-01
IF(NOXCAL) THEN
+ EX(NOSR)=1.981E+00;EX(CRKT)=1.404E+03;EX(XO2 )=1.326E-01
ENDIF
IF(THRAD) THEN
+ EX(SRAD)=3.710E+04
ENDIF