THICKEN_OBSTRUCTIONS = .FALSE.      DROP_VERTICAL_VELOCITY   = 0.20        DROP_HORIZONTAL_VELOCITY = 0.50        DEBUG                = .FALSE.      TIMING               = .FALSE.      FLUSH_FILE_BUFFERS   = .TRUE.C      MISC_LOOP: DO       CALL CHECKREAD('MISC',LU5,IOS) ; IF (IOS.EQ.1) EXIT MISC_LOOP      READ(LU5,MISC,END=23,ERR=24,IOSTAT=IOS)   24 IF (IOS.GT.0) CALL SHUTDOWN('ERROR: Problem with MISC line')      ENDDO MISC_LOOP   23 REWIND(LU5)CC Look to see if DATABASE existsC      IF (DATABASE.EQ.'null' .AND. DATABASE_DIRECTORY.NE.'null') THEN      DATABASE = TRIM(DATABASE_DIRECTORY)//'database4.data'      ENDIFC      NFILES = 1      IF (DATABASE.NE.'null') THEN      INQUIRE(FILE=DATABASE,EXIST=EX)      IF (.NOT.EX) THEN         WRITE(MESSAGE,'(A)') 'ERROR: DATABASE '//TRIM(DATABASE)//     .                        ' not found'         CALL SHUTDOWN(MESSAGE)         ENDIF      NFILES = 2      OPEN(LU80,FILE=DATABASE,FORM='FORMATTED',STATUS='OLD',ERR=190,     .     ACTION='READ',IOSTAT=IOS)  190 IF (IOS.GT.0) CALL SHUTDOWN('ERROR: Problem opening DATABASE')      ENDIFC      H0    = 0.5*(U0**2+V0**2+W0**2)      TMPA  = TMPA + TMPM      TMPA4 = TMPA**4C      IF (TMPO.GT.-1000.) THEN         TMPO = TMPO + TMPM      ELSE          TMPO = TMPA      ENDIFCC Humidity (40% by default, but limited for high temps)C      IF (HUMIDITY.LT.0.)     .HUMIDITY=0.4*MIN(1._EB,EXP(-4890.5*(1/TMPA-1/373.15))/0.0812)CC MiscellaneousC      IF (DENSITY.GT.0.) MW = TMPA*DENSITY*R0/PINF      MW_BACKGROUND = MW      MU_USER(0) = VISCOSITY      K_USER(0)  = THERMAL_CONDUCTIVITY      HCH    = C_HORIZONTAL      HCV    = C_VERTICAL      SURFNAME(0) = SURF_DEFAULT      IF (ISOTHERMAL) RADIATION = .FALSE.      C_FORCED = C_FORCED*(1012.)*(1.8E-5)**0.2 / (0.7)**(2./3.)      TEX_ORI = TEXTURE_ORIGIN      IF (GAUGE_TEMPERATURE.LT.0.) THEN         TMP_GAUGE = TMPA      ELSE         TMP_GAUGE = GAUGE_TEMPERATURE + TMPM      ENDIFC      IF (DTSAM_PART.LT.0.) THEN      WPAR = TWFIN/REAL(NFRAMES,EB)      ELSE      WPAR = DTSAM_PART      ENDIFCC Background Pressure StuffC      P0_MAX = PINF + MAX_OVER_PRESSURE*PINFCC Gravity rampC      IRAMPG = 0      NRAMP  = 0      IF (RAMP_G.NE.'null') THEN      NRAMP = NRAMP + 1      RAMPID(NRAMP) = RAMP_G      RAMPTYPE(NRAMP) = 'TIME'      IRAMPG = NRAMP      ENDIFCC Prandtl and Schmidt numbersC      IF (DNS) THEN         BAROCLINIC = .TRUE.         LES = .FALSE.         IF (PR.LT.0.) PR = 0.7         IF (SC.LT.0.) SC = 1.0         ELSE         IF (PR.LT.0.) PR = 0.5         IF (SC.LT.0.) SC = 0.5         ENDIFC      RSC = 1./SC      RPR = 1./PRCC Check for a restart fileC      APPEND = .FALSE.      IF (RESTART) APPEND = .TRUE.CC Min and Max values of species and temperatureC      TMPMIN = MIN(TMPA,TMPO)      TMPMAX = 3000.      YYMIN  = 0.      YYMAX  = 1.C      END SUBROUTINE READ_MISCCC      SUBROUTINE READ_SPECC      REAL(EB) EPSIJ,SIGMA2,AMW,MASS_FRACTION_0,MW,OMEGA,TSTAR,     .         EPSK(0:5),SIG(0:5),SIGMALJ,EPSILONKLJ,NU,MU_N,K_N,D_N0,     .         ZZ,DZZ,TE,WGT,EPSK_N,SIG_N,MW_N      REAL(EB) CP_N2,CP_F,CP_CO2,CP_O2,CP_H2O,XVAP,MAXIMUM_VISIBILITY      REAL(EB) RCON_MF      INTEGER IPC,NN      LOGICAL READ_MIX      INTEGER  ITMP,IYY,NSPEC_READ      CHARACTER(26) FUEL,STATE_SPECIES(10)      NAMELIST /SURF/ VBC,TMPIGN,TMPWAL,TMPWAL0,ALPHA,C_P,KS,DELTA,     .                C_DELTA_RHO,MASS_FRACTION,VEL,VEL_T,NPPC,     .                E_COEFFICIENT,HEAT_FLUX,WALL_POINTS,     .                TIGN,TAU_Q,TAU_V,RAMP_Q,SURFACE_DENSITY,TAU_MF,     .                RAMP_MF,PART_ID,RAMP_V,VOLUME_FLUX,RAMP_KS,     .                RAMP_KS_CHAR,RAMP_C_P_CHAR,     .                PROFILE,PLE,Z0,ID,AKA,MASS_FLUX,RAMP_C_P,     .                FYI,POROSITY,PAPER_MODEL,     .                BACKING,TMP_BACK,HRRPUA,EMISSIVITY,PHASE,RADIUS,     .                HEAT_OF_VAPORIZATION,DENSITY,IN_DEPTH_COEFFICIENT,     .                HEAT_OF_ABLATION,ABLATION_TEMPERATURE,     .                ABLATION_RATE,ADIABATIC,HEAT_OF_COMBUSTION,     .                BURNING_RATE_MAX,PARTICLES,PARTICLE_COLOR,     .                TEXTURE_MAP,TEXTURE_WIDTH,TEXTURE_HEIGHT,RGB,RGB4,     .                MASS_FLUX_CRITICAL,A,     .                E,MOISTURE_FRACTION,FUEL_FRACTION,     .                CHAR_DENSITY, C_P_CHAR, KS_CHAR,BURN_AWAY,LEAKING,     .                DX_SOLID,EXTERNAL_FLUX,TMPEVAP      NAMELIST /SPEC/ MASS_FRACTION_0,MW,FYI,ID,SIGMALJ,     .                EPSILONKLJ,NU,DENSITY,THERMAL_CONDUCTIVITY,     .                VISCOSITY,DIFFUSION_COEFFICIENT      NAMELIST /REAC/ E,BOF,XNO,XNF,DELTAH,FYI,FUEL,EPUMO2,ID,     .                RADIATIVE_FRACTION,HRRPUA_SHEET,TMP_LOWER,     .                Y_O2_INFTY,Y_F_INLET,NU_FUEL,NU_O2,     .                NU_H2O,NU_CO2,MW_FUEL,FUEL_N2,     .                H2_YIELD,SOOT_YIELD,CO_YIELD,     .                MASS_EXTINCTION_COEFFICIENT,DTSAM,     .                X_O2_LL,Z_CONSTANT,VISIBILITY_FACTOR,     .                MAXIMUM_VISIBILITY,     .                CRITICAL_FLAME_TEMPERATURECC Scan all active SURF lines to see if mixture fraction combustionC model is desiredC      IWATER  = 0      IFUEL   = 0      IOXYGEN = 0           ICO2    = 0      MIXTURE_FRACTION = .FALSE.C      IF (FUEL_EVAPORATION) MIXTURE_FRACTION = .TRUE.C      N = -1      SURF_LOOP: DO       N = N+1      IF (N.GT.NBT) EXIT SURF_LOOPC      AKA_NAME(N) = 'null'      IF (SURFNAME(N).EQ.'INERT') CYCLE SURF_LOOPC      FILE_LOOP: DO ITER=1,NFILESC      IF (ITER.EQ.1) LUDUM = LU5      IF (ITER.EQ.2) LUDUM = LU80      REWIND(LUDUM)C      LOOP1: DO      CALL CHECKREAD('SURF',LUDUM,IOS) ; IF (IOS.EQ.1) CYCLE FILE_LOOP      CALL SET_SURF_DEFAULTS      READ(LUDUM,SURF,ERR=34,IOSTAT=IOS)C      NN = -1      FIND_BACKING: DO      NN = NN+1      IF (NN.GT.NBT) EXIT FIND_BACKING      IF (BACKING.EQ.SURFNAME(NN)) EXIT FIND_BACKING      IF (NN.EQ.NBT .AND. BACKING.NE.'VOID' .AND. BACKING.NE.'EXPOSED'     .    .AND. BACKING.NE.'INSULATED') THEN         NBT = NBT+1         SURFNAME(NBT) = BACKING         EXIT FIND_BACKING         ENDIF      ENDDO FIND_BACKING      IF (ID.EQ.SURFNAME(N) .OR. ID.EQ.AKA_NAME(N)) THEN         IF (AKA.NE.'null') THEN             AKA_NAME(N) = AKA             CYCLE LOOP1             ENDIF         IF (HRRPUA.GT.0. .OR. HEAT_OF_VAPORIZATION.GT.0.) THEN             MIXTURE_FRACTION = .TRUE.             ENDIF          EXIT FILE_LOOP         ENDIF   34 IF (IOS.GT.0) THEN         WRITE(MESSAGE,'(A)') 'ERROR: Problem with SURF '//     .       TRIM(SURFNAME(N))//' or some SURF line preceding it'         CALL SHUTDOWN(MESSAGE)         ENDIF      ENDDO LOOP1C      ENDDO FILE_LOOPC      IF (ID.NE.SURFNAME(N) .AND. ID.NE.AKA_NAME(N)) THEN      WRITE(MESSAGE,'(A)') 'ERROR: SURF '//TRIM(SURFNAME(N))//     .                     ' not found'      CALL SHUTDOWN(MESSAGE)      ENDIFC      REWIND(LU5)      IF (DATABASE.NE.'null') REWIND(LU80)C      ENDDO SURF_LOOPCC Check SPEC lines for errors and determine number of extra speciesC      ALLOCATE(SPECIES_ID(0:5),STAT=IZERO)      CALL ChkMemErr('READ','SPECIES_ID',IZERO)      SPECIES_ID(0) = BACKGROUND_SPECIES      SPECIES_ID(1:5) = 'nullspecies'C      NSPEC = 0      READ_MIX = .FALSE.C      SPECLOOP: DO      CALL CHECKREAD('SPEC',LU5,IOS) ; IF (IOS.EQ.1) EXIT SPECLOOP      ID = 'nullspecies'      NU = 0.      READ(LU5,NML=SPEC,END=29,ERR=30,IOSTAT=IOS)C      IF (ID.EQ.'nullspecies') THEN         WRITE(MESSAGE,'(A,I2,A)') 'ERROR: Species',N,     .                             'needs a name (ID=...)'         CALL SHUTDOWN(MESSAGE)         ENDIFC      NSPEC = NSPEC + 1      SPECIES_ID(NSPEC) = IDC      IF (NU.NE.0.) MIXTURE_FRACTION = .FALSE.      IF (ID.EQ.'MIXTURE_FRACTION') THEN         READ_MIX         = .TRUE.         MIXTURE_FRACTION = .TRUE.         IFUEL            = NSPEC         ENDIF      IF (ID.EQ.'WATER VAPOR') IWATER = NSPEC      IF (ID.EQ.'OXYGEN') IOXYGEN = NSPEC      IF (ID.EQ.'CARBON DIOXIDE') ICO2 = NSPECC   30 IF (IOS.GT.0) THEN         WRITE(MESSAGE,'(A,I2)')      .       'ERROR: Problem with SPECies number',NSPEC+1         CALL SHUTDOWN(MESSAGE)         ENDIF      ENDDO SPECLOOP   29 REWIND(LU5)C      NSPEC_READ = NSPECCC Check if MIXTURE_FRACTION needs to be added to SPEC listC      IF (MIXTURE_FRACTION .AND. IFUEL.EQ.0) THEN         NSPEC = NSPEC + 1         SPECIES_ID(NSPEC) = 'MIXTURE_FRACTION'         IFUEL = NSPEC         ENDIFCC Check if water vapor should be included in the calculationC      IF (WATER_EVAPORATION .AND. IWATER.EQ.0) THEN         NSPEC = NSPEC + 1         SPECIES_ID(NSPEC) = 'WATER VAPOR'         IWATER = NSPEC         ENDIFCC Allocate species-related arraysC      ALLOCATE(NUN(0:NSPEC),STAT=IZERO)      CALL ChkMemErr('READ','NUN',IZERO)C      ALLOCATE(YY0(0:NSPEC),STAT=IZERO)      CALL ChkMemErr('READ','YY0',IZERO)      ALLOCATE(RCON(0:NSPEC),STAT=IZERO)      CALL ChkMemErr('READ','RCON',IZERO)C      IF (MIXTURE_FRACTION) THEN      ALLOCATE(CP(0:100,500),STAT=IZERO)      CALL ChkMemErr('READ','CP',IZERO)      ALLOCATE(RCP(0:100,500),STAT=IZERO)      CALL ChkMemErr('READ','RCP',IZERO)      ALLOCATE(HH(0:100,0:500),STAT=IZERO)      CALL ChkMemErr('READ','HH',IZERO)      ALLOCATE(RCON_STATE(0:10),STAT=IZERO)      CALL ChkMemErr('READ','RCON_STATE',IZERO)      ALLOCATE(MWN(0:10),STAT=IZERO)      CALL ChkMemErr('READ','MWN',IZERO)      ELSE      ALLOCATE(CP(0:NSPEC,500),STAT=IZERO)      CALL ChkMemErr('READ','CP',IZERO)      ALLOCATE(RCP(0:NSPEC,500),STAT=IZERO)      CALL ChkMemErr('READ','RCP',IZERO)      ALLOCATE(HH(0:NSPEC,0:500),STAT=IZERO)      CALL ChkMemErr('READ','HH',IZERO)      ALLOCATE(MWN(0:NSPEC),STAT=IZERO)      CALL ChkMemErr('READ','MWN',IZERO)      ENDIFC      YY0     = 0.      MWN     = 0.      MWN(0)  = MW_BACKGROUND      EPSK    = 0.      SIG     = 0.      NUN     = 0.CC Initialize variables for MIX_FRAC and WATER_VAPOR, if necessaryC      IF (MIXTURE_FRACTION .AND. .NOT.READ_MIX) THEN         YY0(IFUEL) = 0.         MWN(IFUEL) = 0.         EPSK(IFUEL) = 0.         SIG(IFUEL)  = 0.         NUN(IFUEL)  = 0.         ENDIFC      IF (IWATER.GT.0) THEN         LP=>LAGRANGIAN(NPC)         XVAP  = MIN(1._EB,EXP(LP%H_V_0*18./R0*(1./LP%TMP_V-1./     .                     MIN(TMPA,LP%TMP_V))))         YY0(IWATER) = HUMIDITY*XVAP/(29./18. + (1.-29./18.)*XVAP)         MWN(IWATER) = 0.         EPSK(IWATER) = 0.         SIG(IWATER)  = 0.         NUN(IWATER)  = 0.         ENDIFCC Read SPEC info for species other than MIXTURE_FRACTION and WATER VAPORC      SPEC_LOOP: DO N=1,NSPEC_READC      MASS_FRACTION_0 = -1.      MW  = 0.      DENSITY = -1.      EPSILONKLJ = 0.      SIGMALJ    = 0.      NU         = 0.      THERMAL_CONDUCTIVITY  = -1.      VISCOSITY     = -1.      DIFFUSION_COEFFICIENT = -1.C      CALL CHECKREAD('SPEC',LU5,IOS) ; IF (IOS.EQ.1) CYCLE SPEC_LOOP      READ(LU5,SPEC,END=31,ERR=32)C      IF (MASS_FRACTION_0.LT.0.) THEN         YY0(N) = 0.          IF (ID.EQ.'OXYGEN') YY0(N) = 0.23         ELSE         YY0(N) = MASS_FRACTION_0         ENDIFC      IF (DENSITY.GT.0.) MW = TMPA*DENSITY*R0/PINF      MU_USER(N) = VISCOSITY      K_USER(N) = THERMAL_CONDUCTIVITY      D_USER(N) = DIFFUSION_COEFFICIENT      MWN(N) = MW      EPSK(N) = EPSILONKLJ      SIG(N)  = SIGMALJ      NUN(N)  = NU   32 ENDDO SPEC_LOOP   31 REWIND(LU5)CC Set default REACtion parametersC      CALL SET_REAC_DEFAULTSCC Search for user specified REACtion parametersC      REAC_LOOP: DO ITER=1,NFILESC      IF (ITER.EQ.2 .AND. REACTION.EQ.'null') EXIT REAC_LOOPC      IF (ITER.EQ.1) LUDUM = LU5      IF (ITER.EQ.2) LUDUM = LU80      REWIND(LUDUM)C      SCANLOOP: DO      CALL CHECKREAD('REAC',LUDUM,IOS) ; IF (IOS.EQ.1) CYCLE REAC_LOOP      CALL SET_REAC_DEFAULTS      READ(LUDUM,REAC,ERR=434,IOSTAT=IOS)      IF (ID.EQ.REACTION) EXIT REAC_LOOP  434 IF (IOS.GT.0) THEN         WRITE(MESSAGE,'(A)') 'ERROR: Problem with REAC '//     .              TRIM(REACTION)//' or some REAC line preceding it'         CALL SHUTDOWN(MESSAGE)         ENDIF      ENDDO SCANLOOPC      ENDDO REAC_LOOPC      IF (ID.NE.REACTION .AND. REACTION.NE.'null') THEN      WRITE(MESSAGE,'(A,A,A)')      .    'ERROR: REAC ',TRIM(REACTION),' not found'      CALL SHUTDOWN(MESSAGE)      ENDIFC      REWIND(LU5)      IF (DATABASE.NE.'null') REWIND(LU80)C