MODULE READ C      USE PREC      USE VARS      USE CONS      USE TRAN      USE PACKERC      IMPLICIT NONEC      PRIVATE      PUBLIC READ_DATA      CHARACTER(60) FYI      CHARACTER(30) QUANTITY,MAKE,LABEL,CB      CHARACTER(30) COLOR_QUANTITY      CHARACTER(80) MESSAGE      CHARACTER(6)  PHASE      CHARACTER(26) PART_ID,SURF_ID,SURF_IDS(3),SURF_ID6(6),     .              ID,AKA,AKA_NAME(0:200),     .              RAMP_MF(0:20),RAMP_Q,RAMP_V,RAMP_G,SURF_DEFAULT,     .              BACKGROUND_SPECIES,REACTION,RAMP_KS,RAMP_C_P,     .              RAMP_KS_CHAR,RAMP_C_P_CHAR,EVAC_SURF_DEFAULT      CHARACTER(20) PROFILE,BACKING,PARTICLE_COLOR      LOGICAL SUCCESS,BNDF_FACE(-3:3),BNDF_BLOCK,EX,     .        THICKEN_OBSTRUCTIONS,PARTICLES      LOGICAL PAPER_MODEL,BNDF_DEFAULT,ADIABATIC,BAD,     .        BURN_AWAY,LEAKING,OUTLINE      REAL(EB) KS,TAU_MF(0:20)      REAL(EB) XYZ(3),XB(6),VEL_T(2),TEXTURE_ORIGIN(3)      CHARACTER(60) TEXTURE_MAP      REAL(EB) TM,TAU_Q,TAU_V,TIGN,F,HRRPUA,EMISSIVITY,TEXTURE_WIDTH,     .         TEXTURE_HEIGHT,     .         E_COEFFICIENT,RADIATIVE_FRACTION,VOLUME_FLUX,     .         C_HORIZONTAL,C_VERTICAL,     .         ALPHA,C_P,TMPWAL,TMPWAL0,DELTA,VEL,C_DELTA_RHO,VBC,     .         TMPIGN,TMPEVAP,VBC0,DTSAM,PBX,PBY,PBZ,     .         T,MASS_FLUX(0:20),MASS_FRACTION(0:20),XI,ETA,ZETA,     .         Z0,PLE,DELTAH,SURFACE_DENSITY,     .         MW_BACKGROUND,MW,DENSITY,DUMMY,     .         THERMAL_CONDUCTIVITY,VISCOSITY,HEAT_FLUX      REAL(EB) DIFFUSION_COEFFICIENT,MU_USER(0:20),K_USER(0:20),     .         D_USER(0:20),POROSITY,     .         HEAT_OF_VAPORIZATION,HEAT_OF_COMBUSTION,     .         HEAT_OF_ABLATION,ABLATION_TEMPERATURE,ABLATION_RATE,     .         BURNING_RATE_MAX,X_O2_LL,IN_DEPTH_COEFFICIENT,     .         GEN_TO_YLD,HUMIDITY,MW_MIN,MW_MAX,     .         ORIENTATION(3),ROTATION,RGB(3),RGB4(4),RADIUS,     .         MASS_FLUX_CRITICAL,A,     .         E,MOISTURE_FRACTION,FUEL_FRACTION,     .         CHAR_DENSITY, C_P_CHAR, KS_CHAR,     .         CRITICAL_FLAME_TEMPERATURE,DX_SOLID,     .         EXTERNAL_FLUX,TMP_BACK      INTEGER NSPRO,NTCO,II,NSFO,NNN,NR,WALL_POINTS,     .        NPPC,IOR,NSPC,ND,NN,     .        I1,I2,J1,J2,K1,K2,N,I,J,K,IZERO,     .        IOS,NVO,LUDUM,ITER,NFILES,NITER      TYPE (MESH_TYPE), POINTER :: M      TYPE(OBSTRUCTION_TYPE), POINTER :: OB,OB2      TYPE (VENTS_TYPE), POINTER :: VT      TYPE(LAGRANGIAN_TYPE), POINTER :: LPCC      CONTAINSCC      SUBROUTINE READ_DATAC      CALL READ_HEAD      CALL READ_GRID      CALL READ_PDIM      CALL READ_TRAN      CALL READ_TIME      CALL READ_MISC      CALL READ_OBST      CALL READ_VENT      CALL READ_PART      CALL READ_TREE      CALL READ_SPEC      CALL READ_SURF      CALL READ_INIT      CALL READ_RAMP      CALL READ_SPRK      CALL READ_HEAT      CALL READ_SMOD      CALL READ_THCP      CALL READ_PL3D      CALL READ_SLCF      CALL READ_ISOF      CALL READ_BNDFC      END SUBROUTINE READ_DATACCC      SUBROUTINE READ_HEADC      NAMELIST /HEAD/ TITLE,CHID,FYIC      CHID    = 'output'      TITLE   = '      'C      HEAD_LOOP: DO      CALL CHECKREAD('HEAD',LU5,IOS) ; IF (IOS.EQ.1) EXIT HEAD_LOOP      READ(LU5,HEAD,END=13,ERR=14,IOSTAT=IOS)   14 IF (IOS.GT.0) CALL SHUTDOWN('ERROR: Problem with HEAD line')      ENDDO HEAD_LOOP   13 REWIND(LU5)C      CLOOP: DO I=1,39      IF (CHID(I:I).EQ.'.')     .   CALL SHUTDOWN('ERROR: No periods allowed in CHID')      IF (CHID(I:I).EQ.' ') EXIT CLOOP      ENDDO CLOOPC      INQUIRE(FILE=TRIM(CHID)//'.stop',EXIST=EX)      IF (EX) THEN         WRITE(MESSAGE,'(A,A,A)') "ERROR: Remove the file, ",     .   TRIM(CHID)//'.stop',", from the current directory"         CALL SHUTDOWN(MESSAGE)         ENDIFC      END SUBROUTINE READ_HEADCC      SUBROUTINE READ_GRIDC      INTEGER :: IBAR,JBAR,KBAR,NM,POISSON_BC(6)      LOGICAL :: EVACUATION, EVAC_HUMANS      NAMELIST /GRID/ IBAR,JBAR,KBAR,FYI,ID,SYNCHRONIZE,     .     EVACUATION,EVAC_HUMANS,POISSON_BC      TYPE (MESH_TYPE), POINTER :: MC      NMESHES = 0C      GRID_LOOP: DO      CALL CHECKREAD('GRID',LU5,IOS) ; IF (IOS.EQ.1) EXIT GRID_LOOP      READ(LU5,GRID,END=15,ERR=16,IOSTAT=IOS)      NMESHES = NMESHES + 1   16 IF (IOS.GT.0) CALL SHUTDOWN('ERROR: Problem with GRID line')      ENDDO GRID_LOOP   15 REWIND(LU5)C      ALLOCATE(MESH(NMESHES),STAT=IZERO)      CALL ChkMemErr('READ','MESH',IZERO)      ALLOCATE(MESH_NAME(NMESHES),STAT=IZERO)      CALL ChkMemErr('READ','MESH_NAME',IZERO)      ALLOCATE(TUSED(N_TIMERS,NMESHES),STAT=IZERO)      CALL ChkMemErr('READ','TUSED',IZERO)      ALLOCATE(SYNC_TIME_STEP(NMESHES),STAT=IZERO)      CALL ChkMemErr('READ','SYNC_TIME_STEP',IZERO)      SYNC_TIME_STEP = .FALSE.      ALLOCATE(EVACUATION_ONLY(NMESHES),STAT=IZERO)      CALL ChkMemErr('READ','EVACUATION_ONLY',IZERO)      EVACUATION_ONLY = .FALSE.      ALLOCATE(EVACUATION_GRID(NMESHES),STAT=IZERO)      CALL ChkMemErr('READ','EVACUATION_GRID',IZERO)      EVACUATION_GRID = .FALSE.      ALLOCATE(PBC(6,NMESHES),STAT=IZERO)      CALL ChkMemErr('READ','PBC',IZERO)CC Read and store GRID dimensionsC      MESH_LOOP: DO NM=1,NMESHES      IBAR=10 ; JBAR=10 ; KBAR=10      TWO_D = .FALSE.      ID = 'null'      SYNCHRONIZE = .FALSE.      EVACUATION  = .FALSE.      EVAC_HUMANS = .FALSE.      POISSON_BC  = -1      WRITE(MESH_NAME(NM),'(A,I3)') 'MESH',NM      CALL CHECKREAD('GRID',LU5,IOS) ; IF (IOS.EQ.1) EXIT MESH_LOOP      READ(LU5,GRID,END=115)      M => MESH(NM)      M%IBAR = IBAR      M%JBAR = JBAR      M%KBAR = KBAR      M%NEWC = 2*IBAR*JBAR+2*IBAR*KBAR+2*JBAR*KBAR      IF (SYNCHRONIZE) SYNC_TIME_STEP(NM)  = .TRUE.      IF (EVACUATION)  EVACUATION_ONLY(NM) = .TRUE.      IF (EVAC_HUMANS) EVACUATION_GRID(NM) = .TRUE.      IF (JBAR.EQ.1) TWO_D = .TRUE.      IF (TWO_D .AND. JBAR.NE.1) THEN         WRITE(MESSAGE,'(A)') 'ERROR: JBAR must be 1 for all grids'         CALL SHUTDOWN(MESSAGE)         ENDIF      IF (ID.NE.'null') MESH_NAME(NM) = ID      PBC(:,NM) = POISSON_BC(:)      ENDDO MESH_LOOP  115 REWIND(LU5)CC Start the timing arraysC      TUSED      = 0.      TUSED(1,:) = SECOND()C      END SUBROUTINE READ_GRIDCC      SUBROUTINE READ_PDIMC      REAL(EB) :: XBAR,YBAR,ZBAR,XBAR0,YBAR0,ZBAR0      INTEGER :: NMC      REAL(EB) :: RBAR0,RBAR      NAMELIST /PDIM/ XBAR,YBAR,ZBAR,XBAR0,YBAR0,ZBAR0,FYI,RBAR0,RBARC      MESH_LOOP: DO NM=1,NMESHESC      RBAR0 = 0. ; RBAR = -1.      XBAR0 = 0. ; XBAR =  1.      YBAR0 = 0. ; YBAR =  1.      ZBAR0 = 0. ; ZBAR =  1.      CYLINDRICAL = .FALSE.C      CALL CHECKREAD('PDIM',LU5,IOS) ; IF (IOS.EQ.1) EXIT MESH_LOOP      READ(LU5,PDIM,END=19,ERR=20,IOSTAT=IOS)    20 IF (IOS.GT.0) CALL SHUTDOWN('ERROR: Problem with PDIM line')C      M => MESH(NM)C      IF (RBAR.GT.0.) THEN         CYLINDRICAL = .TRUE.         XBAR        = RBAR         XBAR0       = RBAR0         ENDIFC      M%XS    = XBAR0      M%XF    = XBAR      M%YS    = YBAR0      M%YF    = YBAR      M%ZS    = ZBAR0      M%ZF    = ZBAR      M%DXI   = (XBAR-XBAR0)/REAL(M%IBAR,EB)      M%DETA  = (YBAR-YBAR0)/REAL(M%JBAR,EB)      M%DZETA = (ZBAR-ZBAR0)/REAL(M%KBAR,EB)      M%RDXI  = 1./M%DXI      M%RDETA = 1./M%DETA      M%RDZETA= 1./M%DZETA      M%IBM1  = M%IBAR-1      M%JBM1  = M%JBAR-1      M%KBM1  = M%KBAR-1      M%IBP1  = M%IBAR+1      M%JBP1  = M%JBAR+1      M%KBP1  = M%KBAR+1C      ENDDO MESH_LOOP   19 REWIND(LU5)C      END SUBROUTINE READ_PDIMCC      SUBROUTINE READ_TRANCC Compute the polynomial transform function for the vertical coordinateC      REAL(EB), ALLOCATABLE, DIMENSION(:,:) :: A,XX      INTEGER, ALLOCATABLE, DIMENSION(:,:) :: ND      REAL(EB)  PC,CC,COEF      INTEGER  IEXP,IC,IDERIV,N,K,IERROR,IOS,I,MESH_NUMBER,     .         NIPX,NIPY,NIPZ,NIPXS,NIPYS,NIPZS,NIPXF,NIPYF,NIPZF,NM      TYPE (MESH_TYPE), POINTER :: M      TYPE (TRAN_TYPE), POINTER :: T      NAMELIST /TRNX/ IDERIV,CC,PC,FYI,MESH_NUMBER      NAMELIST /TRNY/ IDERIV,CC,PC,FYI,MESH_NUMBER      NAMELIST /TRNZ/ IDERIV,CC,PC,FYI,MESH_NUMBERCC Scan the input file, counting the number of NAMELIST entriesC      ALLOCATE(TRANS(NMESHES))C      MESH_LOOP: DO NM=1,NMESHES      M => MESH(NM)      T => TRANS(NM)C      DO N=1,3      T%NOC(N) = 0      TRNLOOP: DO      IF (N.EQ.1) THEN         CALL CHECKREAD('TRNX',LU5,IOS) ; IF (IOS.EQ.1) EXIT TRNLOOP         MESH_NUMBER = 1         READ(LU5,NML=TRNX,END=17,ERR=18,IOSTAT=IOS)         IF (MESH_NUMBER.NE.NM) CYCLE TRNLOOP         ENDIF      IF (N.EQ.2) THEN         CALL CHECKREAD('TRNY',LU5,IOS) ; IF (IOS.EQ.1) EXIT TRNLOOP         MESH_NUMBER = 1         READ(LU5,NML=TRNY,END=17,ERR=18,IOSTAT=IOS)         IF (MESH_NUMBER.NE.NM) CYCLE TRNLOOP         ENDIF      IF (N.EQ.3) THEN         CALL CHECKREAD('TRNZ',LU5,IOS) ; IF (IOS.EQ.1) EXIT TRNLOOP         MESH_NUMBER = 1         READ(LU5,NML=TRNZ,END=17,ERR=18,IOSTAT=IOS)         IF (MESH_NUMBER.NE.NM) CYCLE TRNLOOP         ENDIF      T%NOC(N) = T%NOC(N) + 1   18 IF (IOS.GT.0) CALL SHUTDOWN('ERROR: Problem with TRN* line')      ENDDO TRNLOOP   17 REWIND(LU5)      ENDDOC      T%NOCMAX = MAX(T%NOC(1),T%NOC(2),T%NOC(3))      ALLOCATE(A(T%NOCMAX+1,T%NOCMAX+1))      ALLOCATE(XX(T%NOCMAX+1,3))      ALLOCATE(ND(T%NOCMAX+1,3))      ALLOCATE(T%C1(0:T%NOCMAX+1,3))      T%C1               = 0.      T%C1(1,1:3)        = 1.      ALLOCATE(T%C2(0:T%NOCMAX+1,3))      ALLOCATE(T%C3(0:T%NOCMAX+1,3))      ALLOCATE(T%CCSTORE(T%NOCMAX,3))      ALLOCATE(T%PCSTORE(T%NOCMAX,3))      ALLOCATE(T%IDERIVSTORE(T%NOCMAX,3))C      T%ITRAN  = 0C      DO IC=1,3      NLOOP:  DO N=1,T%NOC(IC)      IDERIV = -1      IF (IC.EQ.1) THEN         LOOP1: DO         CALL CHECKREAD('TRNX',LU5,IOS) ; IF (IOS.EQ.1) EXIT NLOOP         MESH_NUMBER = 1         READ(LU5,TRNX,END=1,ERR=2)         IF (MESH_NUMBER.EQ.NM) EXIT LOOP1         ENDDO LOOP1         ENDIF      IF (IC.EQ.2) THEN         LOOP2: DO         CALL CHECKREAD('TRNY',LU5,IOS) ; IF (IOS.EQ.1) EXIT NLOOP         MESH_NUMBER = 1         READ(LU5,TRNY,END=1,ERR=2)         IF (MESH_NUMBER.EQ.NM) EXIT LOOP2         ENDDO LOOP2         ENDIF      IF (IC.EQ.3) THEN         LOOP3: DO         CALL CHECKREAD('TRNZ',LU5,IOS) ; IF (IOS.EQ.1) EXIT NLOOP         MESH_NUMBER = 1         READ(LU5,TRNZ,END=1,ERR=2)         IF (MESH_NUMBER.EQ.NM) EXIT LOOP3         ENDDO LOOP3         ENDIF      T%CCSTORE(N,IC) = CC      T%PCSTORE(N,IC) = PC      T%IDERIVSTORE(N,IC) = IDERIV      IF (IDERIV.GE.0) T%ITRAN(IC) = 1      IF (IDERIV.LT.0) T%ITRAN(IC) = 2    2 ENDDO NLOOP    1 REWIND(LU5)      ENDDO C      ICLOOP: DO IC=1,3C      IF (T%ITRAN(IC).EQ.1) THENCC If ITRAN=1, do a polynomial transformation C      ND(1,IC)  = 0      SELECT CASE(IC)      CASE(1)         XX(1,IC)    = M%XF-M%XS         T%C1(1,IC)  = M%XF-M%XS      CASE(2)         XX(1,IC)    = M%YF-M%YS         T%C1(1,IC)  = M%YF-M%YS      CASE(3)         XX(1,IC)    = M%ZF-M%ZS         T%C1(1,IC)  = M%ZF-M%ZS      END SELECT       NNLOOP:  DO N=2,T%NOC(IC)+1      IDERIV = T%IDERIVSTORE(N-1,IC)      IF (IC.EQ.1) CC = T%CCSTORE(N-1,IC)-M%XS      IF (IC.EQ.2) CC = T%CCSTORE(N-1,IC)-M%YS      IF (IC.EQ.3) CC = T%CCSTORE(N-1,IC)-M%ZS      IF (IC.EQ.1 .AND. IDERIV.EQ.0) PC = T%PCSTORE(N-1,IC)-M%XS      IF (IC.EQ.2 .AND. IDERIV.EQ.0) PC = T%PCSTORE(N-1,IC)-M%YS      IF (IC.EQ.3 .AND. IDERIV.EQ.0) PC = T%PCSTORE(N-1,IC)-M%ZS      IF (IC.EQ.1 .AND. IDERIV.GT.0) PC = T%PCSTORE(N-1,IC)      IF (IC.EQ.2 .AND. IDERIV.GT.0) PC = T%PCSTORE(N-1,IC)      IF (IC.EQ.3 .AND. IDERIV.GT.0) PC = T%PCSTORE(N-1,IC)        ND(N,IC) = IDERIV        XX(N,IC) = CC      T%C1(N,IC) = PC      ENDDO NNLOOPC      DO K=1,T%NOC(IC)+1      DO N=1,T%NOC(IC)+1      COEF = IFAC(K,ND(N,IC))      IEXP = K-ND(N,IC)      IF (IEXP.LT.0) A(N,K) = 0.      IF (IEXP.EQ.0) A(N,K) = COEF      IF (IEXP.GT.0) A(N,K) = COEF*XX(N,IC)**IEXP      ENDDO      ENDDOC      IERROR = 0      CALL GAUSSJ(A,T%NOC(IC)+1,T%NOCMAX+1,T%C1(1:T%NOCMAX+1,IC),     .            1,1,IERROR)      IF (IERROR.NE.0)     .   CALL SHUTDOWN('ERROR: Problem with grid transformation')C      ENDIFC      IF (T%ITRAN(IC).EQ.2) THENCC If ITRAN=2, do a linear transformationC      T%C1(0,IC) = 0.      T%C2(0,IC) = 0.      DO N=1,T%NOC(IC)      IF (IC.EQ.1) CC = T%CCSTORE(N,IC)-M%XS      IF (IC.EQ.2) CC = T%CCSTORE(N,IC)-M%YS      IF (IC.EQ.3) CC = T%CCSTORE(N,IC)-M%ZS      IF (IC.EQ.1) PC = T%PCSTORE(N,IC)-M%XS      IF (IC.EQ.2) PC = T%PCSTORE(N,IC)-M%YS      IF (IC.EQ.3) PC = T%PCSTORE(N,IC)-M%ZS      T%C1(N,IC) = CC      T%C2(N,IC) = PC      ENDDOC      SELECT CASE(IC)      CASE(1)         T%C1(T%NOC(1)+1,1) = M%XF-M%XS