TEXT(CH4 Statoil Furnace TITLE DISPLAY The case considered is the turbulent combustion of CH4 with O2 in an industrial furnace. A swirling central jet of 'oxidiser'enters the furnace at 552K with the composition of 80% O2 and 20% H2O, and mixes and burns with an annular 'fuel' jet of 25% CH4, 12%CO, 40% H2O & 23% CO2 at a temperature of 1055 K. 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 kinetic calculation is peformed as a post-processing run with a frozen hydrodynamic and temperature fields. ENDDIS IRUNN = 1 ;LIBREF = 0 CHAR(ANS2,CMOD);BOOLEAN(THRAD,KINET,TWOPNT) REAL(PI,RV,RQ,RF,Z1,Z2,Z3,NCZ1,NCZ2,NCZ3,NCY1,NCY2,NCZALL) REAL(ALPHA,SCALE,KEIN,KESWIRL,EPIN,EPSWIRL) REAL(MAXV,MINCELL,RELX); RELX=1.0E-3 PI = 3.1415 * Radius vortex generator RV = 0.095 * Radius quarl RQ = 0.190 * Radius furnace RF = 0.220 * Axial distance from vortex generator to quarl Z1 = 0.142 * Axial distance over quarl Z2 = 0.261 * Axial distance over furnace Z3 = 2.097 * Number of cells in radial and axial dir. NCY1=9;NCY2=2 NCZ1=4;NCZ2=8;NCZ3=15 NCZALL=NCZ1+NCZ2+NCZ3 * U0 er aksiell hastighet ved innloepet * Re = U0*RV/ENUL = 4.7*0.095/1.78E-06 = 250 000 U0 = 4.7 * WMAX er tagentiell hastighet ved innloept, r=Rmax * WMAX = S0*2*U0 hvor S0=0.7 gir WMAX = 6.6 WMAX = 6.6 ALPHA = 55.0*3.1416/180.0 REAL(WIN,WSWIRL,USWIRL) WIN = 7.0;WSWIRL = 3.0;USWIRL = WSWIRL*TAN(ALPHA) KEIN = 5*(WIN**2)/100 EPIN = 0.09**0.75*KEIN**1.5/0.007*RV/2 KESWIRL = 5*(WSWIRL**2)/100 EPSWIRL = 0.09**0.75*KESWIRL**1.5/0.007*RV/2 *MIB REAL(GASCON,TFUEL,TOX,RHOOX,RHOFU,TWAL) REAL(YCH4F,YH2F,YO2F,YH2OF,YCOF,YCO2F,YN2F) REAL(YCH4O,YH2O,YO2O,YH2OO,YCOO,YCO2O,YN2O) GASCON=8314.0;PRESS0=38.E5;TWAL=800. ** Fuel and oxidiser stream inlet values TOX=552;TFUEL=1055 YCH4F=0.25;YH2F=0.0;YO2F=0.0;YH2OF=0.4;YCOF=0.12 YCO2F=0.23;YN2F=0. YCH4O=0.0;YH2O=0.0;YO2O=0.8;YH2OO=0.2;YCOO=0.0 YCO2O=0.0;YN2O=0. MESG( Do you want a KINETics run? (default=N) READVDU(ANS,CHAR,N) IF(:ANS:.EQ.Y) THEN + KINET=T + MESG( Kinetics run selected. ELSE + KINET=F + MESG( Main combustion run selected. ENDIF Group 1. Run Title Group 2. Transience STEADY = T Groups 3, 4, 5 Grid Information * Overall number of cells, RSET(M,NX,NY,NZ,tolerance) RSET(M,1,NCY1+NCY2,NCZ1+NCZ2+NCZ3) * Set overall domain extent: * xulast yvlast zwlast name XSI= 1.000E+00;YSI= RF;ZSI= Z1+Z2+Z3; RSET(D,CHAM,XSI,YSI,ZSI) Group 6. Body-Fitted coordinates BFC=T * Set points XPO= 0.0000E+00;YPO= 0.0000E+00;ZPO= 0.0000E+00; GSET(P,P1,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= RV ;ZPO= 0.0000E+00; GSET(P,P2,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= RF ;ZPO= 0.0000E+00; GSET(P,P3,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= 0.0000E+00;ZPO= Z1 ; GSET(P,P4,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= RV ;ZPO= Z1 ; GSET(P,P5,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= RF ;ZPO= Z1 ; GSET(P,P6,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= 0.0000E+00;ZPO= Z1+Z2 ; GSET(P,P7,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= RQ ;ZPO= Z1+Z2 ; GSET(P,P8,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= RF ;ZPO= Z1+Z2 ; GSET(P,P9,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= 0.0000E+00;ZPO= Z1+Z2+Z3 ; GSET(P,P10,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= RQ ;ZPO= Z1+Z2+Z3 ; GSET(P,P11,XPO,YPO,ZPO) XPO= 0.0000E+00;YPO= RF ;ZPO= Z1+Z2+Z3 ; GSET(P,P12,XPO,YPO,ZPO) * Set lines/arcs GSET(L,L1,P1,P2,NCY1,1.0) GSET(L,L2,P2,P3,NCY2,1.0) GSET(L,L3,P4,P5,NCY1,1.0) GSET(L,L4,P5,P6,NCY2,1.0) GSET(L,L5,P7,P8,NCY1,1.0) GSET(L,L6,P8,P9,NCY2,1.0) GSET(L,L7,P10,P11,NCY1,1.0) GSET(L,L8,P11,P12,NCY2,1.0) GSET(L,L9,P1,P4,NCZ1,1.0) GSET(L,L10,P4,P7,NCZ2,1.0) GSET(L,L11,P7,P10,NCZ3,1.2) GSET(L,L12,P2,P5,NCZ1,1.0) GSET(L,L13,P5,P8,NCZ2,1.0) GSET(L,L14,P8,P11,NCZ3,1.2) GSET(L,L15,P3,P6,NCZ1,1.0) GSET(L,L16,P6,P9,NCZ2,1.0) GSET(L,L17,P9,P12,NCZ3,1.2) * Set frames GSET(F,F1,P1,-,P2,-,P5,-,P4,-) GSET(F,F2,P2,-,P3,-,P6,-,P5,-) GSET(F,F3,P4,-,P5,-,P8,-,P7,-) GSET(F,F4,P5,-,P6,-,P9,-,P8,-) GSET(F,F5,P7,-,P8,-,P11,-,P10,-) GSET(F,F6,P8,-,P9,-,P12,-,P11,-) * Match a grid mesh GSET(M,F1,+J+K,1,1,1,TRANS) GSET(M,F2,+J+K,1,NCY1+1,1,TRANS) GSET(M,F3,+J+K,1,1,NCZ1+1,TRANS) GSET(M,F4,+J+K,1,NCY1+1,NCZ1+1,TRANS) GSET(M,F5,+J+K,1,1,NCZ1+NCZ2+1,TRANS) GSET(M,F6,+J+K,1,NCY1+1,NCZ1+NCZ2+1,TRANS) * Copy/Transfer/Block grid planes * rotasjon GSET(C,I2,F,I1,1,NCY1+NCY2,1,NCZALL,RZ,-1.0000E-01,0,0,INC,1) * translasjon GSET(C,I4,F,I1,1,NCY1+NCY2,1,NCZALL,+,1.0000E-01,0,0,INC,1) NONORT = T Group 7. Variables: STOREd,SOLVEd,NAMEd IF(KINET) THEN + STORE(P1,U1,V1,W1,KE,EP) ELSE + SOLVE(P1,U1,V1,W1);SOLUTN(P1,P,P,Y,P,P,P) + TURMOD(KECHEN) ENDIF STORE(ENUT,LEN1);STORE(H1,LEN1,RHO1,TMP1) Group 9. Properties ENUL = 1.780E-06 IF(KINET) THEN + CHEMKIN(SPECIES,CH4,O2,H2,CO,H2O,CO2,CH3,CH2,CH,CH2O,HCO,H,O,OH,H$ O2,H2O2,N2) + STORE(F,SPH1,MMWT,SRAD) ELSE *** START OF EXTENDED SCRS MODEL SETTINGS CFUEL=CH4 COXID=O2 CFP1=CO2 CFP2=H2O NRSYS=-2 INTEGER(NSPEC,NELEM);NSPEC=7;NELEM=4 INTEGER(ICOMB,NCSTEP,NCREAC) MESG( Enter required combustion model MESG( Default: 1 step finite rate EBU MESG( The options are: MESG( EBU1 - 1 step finite-rate EBU MESG( EBU2 - 2 step 2 reaction finite-rate EBU MESG( EBU3 - 2 step 3 reaction finite-rate EBU MESG( BURN - Infinite-rate model MESG( DDEL - Infinite-rate model with Double-Delta PDF READVDU(CMOD,CHAR,EBU1) CASE :CMOD: OF WHEN EBU1,4 + MESG(1 step finite-rate EBU model + MESG(2CH4 + 4O2 > 2CO2 + 4H2O + NCSTEP=1;NCREAC=1;ICOMB=1 + SCRS(SYSTEM,NCSTEP,NCREAC,NELEM,FRATE*) + SCRS(SPECIES,CH4,O2,H2,CO,H2O,CO2,N2) + SCRS(PROP,CHEMKIN,SCH4) + STORE(P1RS) WHEN EBU2,4 + MESG(2 step 2 reactions finite-rate EBU model + MESG( 2CH4 + 3O2 > 2CO + 4H2O + MESG( 2CO + O2 > 2CO2 + NCSTEP=2;NCREAC=2;ICOMB=2 ** remove * from FRATE for weak linearisation + SCRS(SYSTEM,NCSTEP,NCREAC,NELEM,FRATE) + SCRS(SPECIES,CH4,O2,H2,CO,H2O,CO2,N2) + SCRS(PROP,CHEMKIN,STWO) + SCRS(EBU,P,2.0);SCRS(EBU,S1,0.00) + STORE(P1RS,S1RS) WHEN EBU3,4 + MESG(2 step 3 reactions finite-rate EBU model + MESG( 2CH4 + O2 > 2CO + 4H2 + MESG( 2CO + O2 > 2CO2 + MESG( 2H2 + O2 > 2H2O + NCSTEP=2;NCREAC=3;ICOMB=3 + SCRS(SYSTEM,NCSTEP,NCREAC,NELEM,FRATE*) + SCRS(SPECIES,CH4,O2,H2,CO,H2O,CO2,N2) + SCRS(PROP,CHEMKIN,SCRS) + SCRS(EBU,P,2.0);SCRS(EBU,S1,0.0);SCRS(EBU,S2,1.0) + STORE(P1RS,S1RS) WHEN BURN,4 + MESG(Infinite rate model + MESG(2CH4 + 4O2 > 2CO2 + 4H2O + NCSTEP=1;NCREAC=1;ICOMB=4 + SCRS(SYSTEM,NCSTEP,NCREAC,NELEM,FASTC) + SCRS(SPECIES,CH4,O2,H2,CO,H2O,CO2,N2) + SCRS(PROP,CHEMKIN,SCH4) WHEN DDEL,4 + MESG(Infinite-rate model with Double-Delta PDF + MESG(2CH4 + 4O2 > 2CO2 + 4H2O + NCSTEP=1;NCREAC=1;ICOMB=5 + SCRS(SYSTEM,NCSTEP,NCREAC,NELEM,FASTC) + SCRS(SPECIES,CH4,O2,H2,CO,H2O,CO2,N2) + SCRS(PROP,CHEMKIN,SCH4) + SCRS(PDF,DDELTA) ENDCASE STORE(MMWT) IF(ICOMB.EQ.3) THEN + YCH4F=0.2;YH2F=0.05 ENDIF ** Define fuel & oxidiser composition & temperature SCRS(FUIN,YCH4F,YO2F,YH2F,YCOF,YH2OF,YCO2F,YN2F,TFUEL) SCRS(OXIN,YCH4O,YO2O,YH2O,YCOO,YH2OO,YCO2O,YN2O,TOX) *** END OF EXTENDED SCRS MODEL SETTINGS ENDIF IF(KINET) THEN + THRAD=F ELSE + MESG( Thermal radiation required ? (default=N) + READVDU(ANS2,CHAR,N) IF(:ANS2:.EQ.Y) THEN + THRAD=T ELSE + THRAD=F ENDIF IF(THRAD) THEN + REAL(ABSORB,SCAT,SIGMA,EMPW,EMISW,EMISG,EMPG) + ABSORB=0.5;SCAT=0.; EMISG=0.07 + SIGMA=5.6697E-8; EMISW=0.4 + 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 ENDIF Group 11.Initialise Var/Porosity Fields CONPOR(QUARL,0.00,VOLUME,#1,#1,#2,#2,#1,#2) CONPOR(LIP,0.00,CELL,1,1,5,5,1,4) IF(THRAD) THEN + REAL(TGUESS);TGUESS=1000.; FIINIT(SRAD)=0.07*SIGMA*TGUESS**4 + FIINIT(H1)=-2.392E+06 ENDIF IF(KINET) THEN + RESTRT(ALL) + DO KK=CHSOB+6,CHSOB+16 + FIINIT(:KK:)=1.E-10 + ENDDO MESG( Do you want the CHEMKIN solver? (default=N) READVDU(ANS,CHAR,N) IF(:ANS:.EQ.Y) THEN + TWOPNT=T;CHEMKIN(REACT,TWOPNT) ** the following statement prevent termination of EARTH run by the TWOPNT solver + SPEDAT(SET,CHEM,BYPASS,L,T) + MESG( CHEMKIN solver selected. ELSE + TWOPNT=F;CHEMKIN(REACT,PHOENICS) + MESG( PHOENICS solver selected. ENDIF + CHSOC=37.5 + CHEMKIN(PROP,STAT,CHSOC) + RHO1=GRND9;TMP1=0.0 TMP1B=13 ELSE + FIINIT(U1)=1.E-5; FIINIT(V1)=1.E-5; FIINIT(W1)=WIN + FIINIT(KE)=KEIN; FIINIT(EP)=EPIN ENDIF Group 13. Boundary & Special Sources * deaktivert fast legeme rotasjon COVAL(SWIRL:JJ:, UC1 , ONLYMS, USWIRL*JJ/10) IF(KINET) THEN DO JJ = 1,4 + INLET(NOCPCKO:JJ:,LOW,1,1,JJ,JJ,1,1,1,LSTEP) + VALUE(NOCPCKO:JJ:,P1,GRND9); VALUE(NOCPCKO:JJ:,W1,WSWIRL) + VALUE(NOCPCKO:JJ:,UCRT,0.0); VALUE(NOCPCKO:JJ:,VCRT,0.0) + VALUE(NOCPCKO:JJ:,TMP1,TOX) + VALUE(NOCPCKO:JJ:,O2,YO2O); VALUE(NOCPCKO:JJ:,H2O,YH2OO) ENDDO + INLET(NOCPCKF,LOW,1,1,6,9,1,1,1,LSTEP) + VALUE(NOCPCKF,P1,GRND9); VALUE(NOCPCKF,W1,WIN) + VALUE(NOCPCKF,TMP1,TFUEL) + VALUE(NOCPCKF,H2O,YH2OF); VALUE(NOCPCKF,CO,YCOF) + VALUE(NOCPCKF,CH4,YCH4F); VALUE(NOCPCKF,CO2,YCO2F) ELSE DO JJ = 1,4 + INLET(SCRSO:JJ:,LOW,1,1,JJ,JJ,1,1,1,LSTEP) + VALUE(SCRSO:JJ:,P1,GRND3); VALUE(SCRSO:JJ:,U1,USWIRL) + VALUE(SCRSO:JJ:,V1,GRND3); VALUE(SCRSO:JJ:,W1,GRND3) + VALUE(SCRSO:JJ:,UCRT,0.0); VALUE(SCRSO:JJ:,VCRT,0.0) + VALUE(SCRSO:JJ:,WCRT,WSWIRL) + VALUE(SCRSO:JJ:,KE,KESWIRL); VALUE(SCRSO:JJ:,EP,EPSWIRL) + VALUE(SCRSO:JJ:,H1,GRND3) + VALUE(SCRSO:JJ:,F,0.0); VALUE(SCRSO:JJ:,CH4,YCH4O) + VALUE(SCRSO:JJ:,CO,YCOO); VALUE(SCRSO:JJ:,H2,YH2O) ENDDO + INLET(SCRSF,LOW,1,1,6,9,1,1,1,LSTEP) + VALUE(SCRSF,P1,GRND3); VALUE(SCRSF,KE,KEIN) + VALUE(SCRSF,EP,EPIN); VALUE(SCRSF,H1,GRND3) + VALUE(SCRSF,F,1.0); VALUE(SCRSF,CH4,YCH4F); VALUE(SCRSF,CO,YCOF) + VALUE(SCRSF,U1,0.0); VALUE(SCRSF,V1,GRND3); VALUE(SCRSF,W1,GRND3) + VALUE(SCRSF,UCRT,0.0); VALUE(SCRSF,VCRT,0.0);VALUE(SCRSF,WCRT,WIN) + VALUE(SCRSF,H2,YH2F) ENDIF ## Coefficient at outlet controls pressure distribution ** Exit Boundary Conditions PATCH(OUT,HIGH,#1,#1,#1,#2,#3,#3,#1,#1) REAL(CMOUT);CMOUT=1.0 IF(KINET) THEN ** convert CMOUT from s/m**3 to g*s/(kg*cm**3) + COVAL(OUT,P1,CMOUT*1.E-3,0.0) + DO KK=CHSOB,CHSOB+16 + VALUE(OUT,:KK:,0.0) + ENDDO IF(TWOPNT) THEN + DO KK=CHSOB,CHSOB+16 + RESREF(:KK:)=-1.E-12 + ENDIT(:KK:)=1.E-6 + ENDDO ENDIF ELSE + VALUE(OUT,F,0.0);COVAL(OUT,P1,CMOUT,0.0) + DO KK=CHSOB,CHSOB+12 + VALUE(OUT,:KK:,0.0) + ENDDO ENDIF PATCH(OUTER,NWALL,1,NX,NY,NY,#3,#NREGZ,1,1) COVAL(OUTER,KE,GRND2,GRND2);COVAL(OUTER,EP,GRND2,GRND2) COVAL(OUTER,U1,GRND2,0);COVAL(OUTER,W1,GRND2,0) PATCH(QUARLN,NWALL,#1,#1,#1,#1,#1,#2,1,1) COVAL(QUARLN,KE,GRND2,GRND2);COVAL(QUARLN,EP,GRND2,GRND2) COVAL(QUARLN,U1,GRND2,0);COVAL(QUARLN,W1,GRND2,0) PATCH(QUARLL,LWALL,#1,#1,#2,#2,#3,#3,1,1) COVAL(QUARLL,KE,GRND2,GRND2);COVAL(QUARLL,EP,GRND2,GRND2) COVAL(QUARLL,U1,GRND2,0);COVAL(QUARLL,V1,GRND2,0) PATCH(I_W_S,SWALL,1,1,5,5,1,4,1,LSTEP) COVAL(I_W_S,KE,GRND2,GRND2);COVAL(I_W_S,EP,GRND2,GRND2) COVAL(I_W_S,U1,GRND2,0);COVAL(I_W_S,V1,GRND2,0) COVAL(I_W_S,W1,GRND2,0) PATCH(S_W_N,NWALL,1,1,5,5,1,4,1,LSTEP) COVAL(S_W_N,KE,GRND2,GRND2);COVAL(S_W_N,EP,GRND2,GRND2) COVAL(S_W_N,U1,GRND2,0);COVAL(S_W_N,V1,GRND2,0) COVAL(S_W_N,W1,GRND2,0) PATCH(L_W_L,LWALL,1,1,5,5,4,4,1,LSTEP) COVAL(L_W_L,KE,GRND2,GRND2);COVAL(L_W_L,EP,GRND2,GRND2) COVAL(L_W_L,U1,GRND2,0);COVAL(L_W_L,V1,GRND2,0) COVAL(L_W_L,W1,GRND2,0) IF(THRAD) THEN + PATCH(NWALL3R,NORTH,1,NX,#NREGY,#NREGY,#3,#NREGZ,#1,#NREGT) + COVAL(NWALL3R,SRAD,EMISW/(2.0-EMISW),EMPW) ENDIF UCONV = T PATCH(SWRC,CELL,1,NX,1,NY,1,NZ,1,1);COVAL(SWRC,U1,FIXFLU,GRND) PATCH(SWRV,VOLUME,1,NX,1,NY,1,NZ,1,1) COVAL(SWRV,V1,FIXFLU,GRND);COVAL(SWRV,W1,FIXFLU,GRND) Group 15. Terminate Sweeps LSWEEP = 2000 Group 17. Relaxation IF(KINET) THEN IF(TWOPNT) THEN + CHEMKIN(RELAX,1.E-4) ELSE + CHEMKIN(RELAX,1.E-2) + SELREF=T ENDIF ** 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 kinetics run. However, on kinetics restart runs you must set LG(94)=F + RELAX(DEN1,LINRLX,1.E-10) MESG( RE-STARTING from SCRS combustion solution? (default=Y) READVDU(ANS,CHAR,Y) IF(:ANS:.EQ.Y) THEN + LG(94)=T + MESG( RESTRT run from SCRS solution ELSE + LG(94)=F; RESTRT(CH3,CH2,CH,CH2O,HCO,H,O,OH,HO2,H2O2,N2) + MESG( RESTRT run from KINETICS solution ENDIF ELSE + RELAX(P1,LINRLX,0.3); RELAX(U1,FALSDT,RELX); RELAX(V1,FALSDT,RELX) + RELAX(W1,FALSDT,RELX); RELAX(KE,FALSDT,RELX);RELAX(EP,FALSDT,RELX) + RELAX(F,LINRLX,0.5); RELAX(RHO1,LINRLX,0.5); RELAX(H1,LINRLX,0.5) IF(THRAD) THEN + RELAX(SRAD,FALSDT,1.0) ENDIF + KELIN=1 IF(ICOMB.LE.3) THEN + RELAX(CH4,FALSDT,RELX);RELAX(CO,FALSDT,RELX);RELAX(H2,FALSDT,RELX) ENDIF IF(ICOMB.EQ.3) THEN + RELAX(RHO1,LINRLX,0.05) ENDIF ENDIF Group 18. Limits Group 19. EARTH Calls To GROUND Station Group 20. Preliminary Printout OUTPUT(SPH1,P,P,Y,P,P,P) Group 22. Monitor Print-Out IXMON=1;IYMON=4;IZMON=24;TSTSWP=-1;ITABL=3 Group 23.Field Print-Out & Plot Control Group 24. Dumps For Restarts ** Activate thermal radiation after solution has developed; following commands required IF(THRAD) THEN IF(KINET) THEN ELSE + MESG( Enter mode of radiation calculation + MESG( Default: Restart from non-radiative SCRS run + MESG( The options are: + MESG( RAD1 - Re-starting from non-radiative SCRS run + MESG( RAD2 - Re-starting from radiative SCRS run + MESG( RAD3 - Starting from scratch + READVDU(CMOD,CHAR,RAD1) CASE :CMOD: OF WHEN RAD1,4 + MESG( RE-STARTING from non-radiative SCRS combustion run + RESTRT(P1,U1,V1,W1,KE,EP,RHO1,H1,CH4,CO2,H2,N2) + RESTRT(H2O,CO2,CO,O2,F,ENUT,UCRT,VCRT,WCRT,TMP1) + FIINIT(SRAD)=0.07*SIGMA*TGUESS**4 WHEN RAD2,4 + MESG( RE-STARTING from radiative SCRS combustion run + RESTRT(P1,U1,V1,W1,KE,EP,RHO1,H1,CH4,CO2,H2,N2) + RESTRT(H2O,CO2,CO,O2,F,ENUT,UCRT,VCRT,WCRT,TMP1,SRAD) WHEN RAD3,4 + MESG( Starting from scratch + FIINIT(SRAD)=0.07*SIGMA*TGUESS**4 ENDCASE ENDIF ENDIF NZPRIN=1;NYPRIN=1;IZPRF=12;IZPRL=16;NPLT=100 IDISPA=10000;CSG1=SW