wall.f90

FDS为火灾动力学模拟软件源代码,该软件为开源项目,代码语言主要为FORTRAN,可在WINDOWS和LINUX下编译运行,详细说明可参考http://fire.nist.gov/fds/官方网址

      ENDIF
      Y_SUM_W = 0._EB
      DO J=1,N_SPECIES
         IF (SPECIES(J)%MODE==GAS_SPECIES) THEN
            Y_SUM_W = Y_SUM_W + YY_W(IW,J)
         ENDIF
      ENDDO            
      CALL GET_MOLECULAR_WEIGHT2(Y_MF_W ,Z_2,YY_W(IW,I_PROG_F),Y_SUM_W,RSUM_W)
      ! Weighting for wall and gas values
      YY_S    = 0.2*Y_MF_W+0.8*Y_MF_G
      RSUM_S  = 0.2*RSUM_W+0.8*RSUM_G
      HVRG    = REACTION(1)%MW_FUEL*ML%H_R(1)/R0
      PPSURF  = MIN(1._EB,YY_S*R0/(REACTION(1)%MW_FUEL*RSUM_S))
      PPCLAUS = MIN(1._EB,EXP(HVRG*(1./ML%TMP_BOIL(1)-1./WC%TMP_S(1))))
      IF (PPSURF/=PPCLAUS) THEN
         IF (MFLUX==0.AND.PPSURF<PPCLAUS) THEN
            MFLUX = 2.E-5*REACTION(1)%MW_FUEL
            TW(IW)=T
         ELSE
            MFLUX = MFLUX*MIN(1.02_EB,MAX(0.98_EB,PPCLAUS/PPSURF))
         ENDIF
      ENDIF
      IF (WC%TMP_S(1)>ML%TMP_BOIL(1)) THEN
         ! Net flux guess for liquid evap
         QNETF = Q_WATER_F + QRADIN(IW) - QRADOUT(IW) + QCONF(IW)
         MFLUX = MAX(MFLUX,1.02*(QNETF - 2.*(ML%TMP_BOIL(1)-WC%TMP_S(1))/DXF/K_S(1))/ML%H_R(1))
      ENDIF
      MFLUX = MIN(MFLUX,THICKNESS*ML%RHO_S/DT_BC)
      ! CYLINDRICAL and SPHERICAL scaling not implemented
      MASSFLUX(IW,I_FUEL) = MASSFLUX(IW,I_FUEL)  + ML%ADJUST_BURN_RATE*ML%NU_FUEL(1)*MFLUX
      ACTUAL_BURN_RATE(IW)= ACTUAL_BURN_RATE(IW) +                     ML%NU_FUEL(1)*MFLUX
      IF (I_WATER>0) MASSFLUX(IW,I_WATER) = MASSFLUX(IW,I_WATER) + ML%NU_WATER(J)*MFLUX
      DX_GRID = DT_BC*MFLUX/ML%RHO_S
      Q_S(1) = Q_S(1) - MFLUX*ML%H_R(1)/DX_S(1)  ! no improvement (in cone test) if used updated RDX 

      IF (POINT_SHRINK) THEN
         X_S_NEW(1:NWP) = MAX(0._EB,X_S_NEW(1:NWP)-DX_GRID)
         IF (DX_GRID > 0._EB) WC%SHRINKING = .TRUE.
      ENDIF
      
      EXIT MATERIAL_LOOP2     ! Can handle only one LIQUID fuel at the time

   ENDDO MATERIAL_LOOP2

   ! Re-generate grid for shrinking wall
   N_LAYER_CELLS_NEW = 0
   RECOMPUTE_GRID: IF (WC%SHRINKING) THEN
      NWP_NEW = 0
      THICKNESS = 0._EB
      RECOMPUTE = .FALSE.
      I = 0
      LAYER_LOOP: DO NL=1,SF%N_LAYERS
         WC%LAYER_THICKNESS(NL) = X_S_NEW(I+WC%N_LAYER_CELLS(NL))-X_S_NEW(I)
         ! Remove very thin layers
         IF (WC%LAYER_THICKNESS(NL) < 1.0E-6_EB) THEN 
            X_S_NEW(I+WC%N_LAYER_CELLS(NL):NWP) = X_S_NEW(I+WC%N_LAYER_CELLS(NL):NWP)-WC%LAYER_THICKNESS(NL)
            WC%LAYER_THICKNESS(NL) = 0._EB
            N_LAYER_CELLS_NEW(NL)  = 0
         ELSE
            CALL GET_N_LAYER_CELLS(SF%MIN_DIFFUSIVITY(NL),WC%LAYER_THICKNESS(NL), &
               SF%STRETCH_FACTOR,SF%CELL_SIZE_FACTOR,N_LAYER_CELLS_NEW(NL),SMALLEST_CELL_SIZE(NL))
            NWP_NEW = NWP_NEW + N_LAYER_CELLS_NEW(NL)
         ENDIF
         IF ( N_LAYER_CELLS_NEW(NL) .NE. WC%N_LAYER_CELLS(NL)) RECOMPUTE = .TRUE.
         THICKNESS = THICKNESS + WC%LAYER_THICKNESS(NL)
         I = I + WC%N_LAYER_CELLS(NL)
      ENDDO LAYER_LOOP
      DO I = 1,NWP
         IF ( (X_S_NEW(I)-X_S_NEW(I-1)) < 0.5 * SMALLEST_CELL_SIZE(LAYER_INDEX(I))) RECOMPUTE = .TRUE.
      ENDDO
      ! Shrinking wall has gone to zero thickness.
      IF (THICKNESS == 0._EB) THEN
         WC%TMP_S(0:NWP+1)    = MAX(TMPMIN,SF%TMP_BACK)
         TMP_F(IW)            = MIN(TMPMAX,MAX(TMPMIN,SF%TMP_BACK))
         TMP_B(IW)            = MIN(TMPMAX,MAX(TMPMIN,SF%TMP_BACK))
         RHOWAL               = 0.5_EB*(RHO(IIG,JJG,KKG)+RHO_W(IW))
         CP_TERM              = MAX(0._EB,-CP_GAMMA*UW(IW)*RHOWAL)
         QCONF(IW)            = 0._EB
         TMP_W(IW)            = ( (RDN(IW)*KW(IW)-0.5_EB*CP_TERM)*TMP_G+CP_TERM*TMP_F(IW)-QCONF(IW) ) & 
                                /(0.5_EB*CP_TERM+RDN(IW)*KW(IW))
         TMP_W(IW)            = MAX(TMPMIN,TMP_W(IW))
         MASSFLUX(IW,I_FUEL)  = 0._EB
         ACTUAL_BURN_RATE(IW) = 0._EB
         IF (I_WATER>0) MASSFLUX(IW,I_WATER)  = 0._EB
         WC%N_LAYER_CELLS     = 0
         WC%BURNAWAY          = .TRUE.
         I_OBST               = OBST_INDEX_W(IW)
         IF (I_OBST > 0) THEN
            IF (OBSTRUCTION(I_OBST)%CONSUMABLE) OBSTRUCTION(I_OBST)%MASS = -1.
         ENDIF
         CYCLE WALL_CELL_LOOP
      ENDIF

      WC%X_S(0:NWP) = X_S_NEW(0:NWP)
      X_S_NEW = 0._EB
      IF (RECOMPUTE) THEN
         CALL GET_WALL_NODE_COORDINATES(NWP_NEW,SF%N_LAYERS,N_LAYER_CELLS_NEW, &
            SMALLEST_CELL_SIZE(1:SF%N_LAYERS),SF%STRETCH_FACTOR,X_S_NEW(0:NWP_NEW))
         CALL GET_WALL_NODE_WEIGHTS(NWP_NEW,SF%N_LAYERS,N_LAYER_CELLS_NEW,THICKNESS,SF%GEOMETRY, &
            X_S_NEW(0:NWP_NEW),DX_S(1:NWP_NEW),RDX_S(0:NWP_NEW+1),RDXN_S(0:NWP_NEW),&
            DX_WGT_S(0:NWP_NEW),DXF,DXB,LAYER_INDEX)           
         ! Interpolate densities and temperature from old grid to new grid
         CALL GET_INTERPOLATION_WEIGHTS(SF%N_LAYERS,NWP,NWP_NEW,WC%N_LAYER_CELLS,N_LAYER_CELLS_NEW, &
                                    WC%X_S(0:NWP),X_S_NEW(0:NWP_NEW),INT_WGT(1:NWP_NEW,1:NWP))      
         CALL INTERPOLATE_WALL_ARRAY(NWP,NWP_NEW,INT_WGT(1:NWP_NEW,1:NWP),WC%TMP_S(1:NWP))
         WC%TMP_S(NWP_NEW+1) = WC%TMP_S(NWP+1)
         CALL INTERPOLATE_WALL_ARRAY(NWP,NWP_NEW,INT_WGT(1:NWP_NEW,1:NWP),Q_S(1:NWP))
         DO N=1,SF%N_MATL
            ML  => MATERIAL(SF%MATL_INDEX(N))
            CALL INTERPOLATE_WALL_ARRAY(NWP,NWP_NEW,INT_WGT(1:NWP_NEW,1:NWP),WC%RHO_S(1:NWP,N))
         ENDDO
         WC%N_LAYER_CELLS = N_LAYER_CELLS_NEW(1:SF%N_LAYERS)
         NWP = NWP_NEW
         WC%X_S(0:NWP) = X_S_NEW(0:NWP)      ! Note: WC%X_S(NWP+1...) are not set to zero.
      ELSE      
         CALL GET_WALL_NODE_WEIGHTS(NWP,SF%N_LAYERS,N_LAYER_CELLS_NEW,THICKNESS,SF%GEOMETRY, &
            WC%X_S(0:NWP),DX_S(1:NWP),RDX_S(0:NWP+1),RDXN_S(0:NWP),DX_WGT_S(0:NWP),DXF,DXB,LAYER_INDEX)
      ENDIF
      
   ENDIF RECOMPUTE_GRID

   ELSEIF (SF%PYROLYSIS_MODEL==PYROLYSIS_SPECIFIED) THEN

      ! Take off energy corresponding to specified burning rate
      Q_S(1) = Q_S(1) - MASSFLUX(IW,I_FUEL)*SF%H_V/DX_S(1)

   ENDIF PYROLYSIS_MATERIAL_IF

   ! Calculate thermal properties 
 
   K_S     = 0._EB
   RHO_S   = 0._EB
   RHOCBAR = 0._EB
   E_WALL(IW) = 0._EB
   POINT_LOOP3: DO I=1,NWP
      VOLSUM = 0._EB
      MATERIAL_LOOP3: DO N=1,SF%N_MATL
         IF (WC%RHO_S(I,N)==0._EB) CYCLE MATERIAL_LOOP3
         ML  => MATERIAL(SF%MATL_INDEX(N))
         VOLSUM = VOLSUM + WC%RHO_S(I,N)/ML%RHO_S
         IF (ML%K_S>0._EB) THEN  
            K_S(I) = K_S(I) + WC%RHO_S(I,N)*ML%K_S/ML%RHO_S
         ELSE
            ITMP   = NINT(WC%TMP_S(I))
            NR     = -NINT(ML%K_S)
            K_S(I) = K_S(I) + WC%RHO_S(I,N)*EVALUATE_RAMP(WC%TMP_S(I),0._EB,NR)/ML%RHO_S
         ENDIF

         IF (ML%C_S>0._EB) THEN
            RHOCBAR(I) = RHOCBAR(I) + WC%RHO_S(I,N)*ML%C_S
         ELSE
            ITMP   = NINT(WC%TMP_S(I))
            NR     = -NINT(ML%C_S)
            RHOCBAR(I) = RHOCBAR(I) + WC%RHO_S(I,N)*EVALUATE_RAMP(WC%TMP_S(I),0._EB,NR)*C_S_ADJUST_UNITS
         ENDIF
         IF (I.EQ.1) E_WALL(IW) = E_WALL(IW) + WC%RHO_S(I,N)*ML%EMISSIVITY/ML%RHO_S
         RHO_S(I)   = RHO_S(I) + WC%RHO_S(I,N)
         
      ENDDO MATERIAL_LOOP3
      
      K_S(I) = K_S(I)/VOLSUM
      IF (I.EQ.1) E_WALL(IW) = E_WALL(IW)/VOLSUM
   ENDDO POINT_LOOP3
  
   K_S(0)     = K_S(1)     ! Calculate average K_S between at grid cell boundaries. Store result in K_S
   K_S(NWP+1) = K_S(NWP)
   DO I=1,NWP-1 
      K_S(I)  = 1.0_EB / ( DX_WGT_S(I)/K_S(I) + (1.-DX_WGT_S(I))/K_S(I+1) )
   ENDDO

   !RCP_W for part
   IF (N_EVAP_INDICIES>0) THEN
      RCP_W(IW) = RHOCBAR(1)*DX_S(1)
   ENDIF


   ! Calculate internal radiation
      
   IF (SF%INTERNAL_RADIATION) THEN

      KAPPA_S = 0._EB
      DO I=1,NWP
         VOLSUM = 0._EB
         DO N=1,SF%N_MATL
            IF (WC%RHO_S(I,N)==0._EB) CYCLE
            ML  => MATERIAL(SF%MATL_INDEX(N))
            VOLSUM = VOLSUM + WC%RHO_S(I,N)/ML%RHO_S
            KAPPA_S(I) = KAPPA_S(I) + WC%RHO_S(I,N)*ML%KAPPA_S/ML%RHO_S
         ENDDO
         KAPPA_S(I) = 2.*KAPPA_S(I)*DX_S(I)/VOLSUM    ! kappa = 2*dx*kappa
      ENDDO
      ! solution inwards
      RFLUX_UP = QRADIN(IW) + (1.-E_WALL(IW))*QRADOUT(IW)
      DO I=1,NWP
         RFLUX_DOWN =  ( RFLUX_UP + KAPPA_S(I)*SIGMA*WC%TMP_S(I)**4 ) / (1. + KAPPA_S(I))
         Q_S(I) = Q_S(I) + (RFLUX_UP - RFLUX_DOWN)/DX_S(I)
         RFLUX_UP = RFLUX_DOWN
      ENDDO
      IF (SF%BACKING==EXPOSED) THEN
         IF (BOUNDARY_TYPE(IWB)==SOLID_BOUNDARY) QRADOUT(IWB) = E_WALLB*RFLUX_UP
      ENDIF
      ! solution outwards   
      RFLUX_UP = QRADINB + (1.-E_WALLB)*RFLUX_UP
      DO I=NWP,1,-1
         RFLUX_DOWN =  ( RFLUX_UP + KAPPA_S(I)*SIGMA*WC%TMP_S(I)**4 ) / (1. + KAPPA_S(I))
         Q_S(I) = Q_S(I) + (RFLUX_UP - RFLUX_DOWN)/DX_S(I)
         RFLUX_UP = RFLUX_DOWN
      ENDDO
      QRADOUT(IW) = E_WALL(IW)*RFLUX_DOWN
   ENDIF

   ! Update the 1-D heat transfer equation 

   DT2_BC = DT_BC
   STEPCOUNT = 1
   ALLOCATE(TMP_W_NEW(0:NWP+1))
   TMP_W_NEW = WC%TMP_S(0:NWP+1)
   WALL_ITERATE: DO 
      ITERATE=.FALSE.
      SUB_TIME: DO N=1,STEPCOUNT
         DXKF   = K_S(0)/DXF
         DXKB   = K_S(NWP)/DXB

         DO I=1,NWP
            BBS(I) = -0.5_EB*DT2_BC*K_S(I-1)*RDXN_S(I-1)*RDX_S(I)/RHOCBAR(I) ! DT_BC->DT2_BC
            AAS(I) = -0.5_EB*DT2_BC*K_S(I)  *RDXN_S(I)  *RDX_S(I)/RHOCBAR(I)
         ENDDO
         DDS(1:NWP) = 1._EB    - AAS(1:NWP) - BBS(1:NWP)
         DO I=1,NWP
            CCS(I) = TMP_W_NEW(I) - AAS(I)*(TMP_W_NEW(I+1)-TMP_W_NEW(I)) + BBS(I)*(TMP_W_NEW(I)-TMP_W_NEW(I-1)) &
                     + DT2_BC*Q_S(I)/RHOCBAR(I)
         ENDDO
         IF (SF%INTERNAL_RADIATION) THEN
            RFACF = 0.25_EB*HEAT_TRANS_COEF(IW)
            RFACB = 0.25_EB*HTCB
         ELSE
            RFACF = 0.25_EB*HEAT_TRANS_COEF(IW)+2._EB*E_WALL(IW)*SIGMA*TMP_F(IW)**3
            RFACB = 0.25_EB*HTCB               +2._EB*E_WALLB*   SIGMA*TMP_B(IW)**3
         ENDIF
         RFACF2 = (DXKF-RFACF)/(DXKF+RFACF)
         RFACB2 = (DXKB-RFACB)/(DXKB+RFACB)
         IF (SF%INTERNAL_RADIATION) THEN
            QDXKF = (HEAT_TRANS_COEF(IW)*(TMP_G   - 0.5_EB*TMP_F(IW)) + Q_WATER_F)/(DXKF+RFACF)
            QDXKB = (HTCB*               (TMP_G_B - 0.5_EB*TMP_B(IW)) + Q_WATER_B)/(DXKB+RFACB)
         ELSE
            QDXKF = (HEAT_TRANS_COEF(IW)*(TMP_G   - 0.5_EB*TMP_F(IW)) + QRADIN(IW) + 3.*E_WALL(IW)*SIGMA*TMP_F(IW)**4 + Q_WATER_F) &
                  /(DXKF+RFACF)
            QDXKB = (HTCB*               (TMP_G_B - 0.5_EB*TMP_B(IW)) + QRADINB    + 3.*E_WALLB   *SIGMA*TMP_B(IW)**4 + Q_WATER_B) &
                  /(DXKB+RFACB)
         ENDIF
         CCS(1)   = CCS(1)   - BBS(1)  *QDXKF
         CCS(NWP) = CCS(NWP) - AAS(NWP)*QDXKB
         DDT(1:NWP) = DDS(1:NWP)
         DDT(1)   = DDT(1)   + BBS(1)  *RFACF2
         DDT(NWP) = DDT(NWP) + AAS(NWP)*RFACB2
         TRIDIAGONAL_SOLVER_1: DO I=2,NWP
            RR     = BBS(I)/DDT(I-1)
            DDT(I) = DDT(I) - RR*AAS(I-1)
            CCS(I) = CCS(I) - RR*CCS(I-1)
         ENDDO TRIDIAGONAL_SOLVER_1
         CCS(NWP)  = CCS(NWP)/DDT(NWP)
         TRIDIAGONAL_SOLVER_2: DO I=NWP-1,1,-1
            CCS(I) = (CCS(I) - AAS(I)*CCS(I+1))/DDT(I)
         ENDDO TRIDIAGONAL_SOLVER_2
         TMP_W_NEW(1:NWP) = MAX(TMPMIN,CCS(1:NWP))
         TMP_W_NEW(0)     = MAX(TMPMIN,TMP_W_NEW(1)  *RFACF2+QDXKF)
         TMP_W_NEW(NWP+1) = MAX(TMPMIN,TMP_W_NEW(NWP)*RFACB2+QDXKB)
         IF (STEPCOUNT==1) THEN
            TOLERANCE = MAXVAL(ABS((TMP_W_NEW-WC%TMP_S(0:NWP+1))/WC%TMP_S(0:NWP+1)), &
               WC%TMP_S(0:NWP+1)>0._EB) ! returns a negative number, if all TMP_S == 0.
            IF (TOLERANCE<0.0_EB) &
            TOLERANCE = MAXVAL(ABS((TMP_W_NEW-WC%TMP_S(0:NWP+1))/TMP_W_NEW), &
               TMP_W_NEW>0._EB) 
            IF (TOLERANCE > 0.2_EB) THEN
               STEPCOUNT = MIN(200,STEPCOUNT * (INT(TOLERANCE/0.2_EB) + 1))
               ITERATE=.TRUE.
               DT2_BC=DT_BC/REAL(STEPCOUNT)
               TMP_W_NEW = WC%TMP_S(0:NWP+1)
            ENDIF
         ENDIF
         TMP_F(IW)       = 0.5_EB*(TMP_W_NEW(0)+TMP_W_NEW(1))
         TMP_F(IW)       = MIN(TMPMAX,MAX(TMPMIN,TMP_F(IW)))
         TMP_B(IW)       = 0.5_EB*(TMP_W_NEW(NWP)+TMP_W_NEW(NWP+1))
         TMP_B(IW)       = MIN(TMPMAX,MAX(TMPMIN,TMP_B(IW)))
      ENDDO SUB_TIME
      IF (.NOT. ITERATE) EXIT WALL_ITERATE
   ENDDO WALL_ITERATE
   WC%TMP_S(0:NWP+1) = TMP_W_NEW
   DEALLOCATE(TMP_W_NEW)

   ! If the surface temperature exceeds the ignition temperature, burn it

   IF (TW(IW) > T ) THEN
      IF (TMP_F(IW) >= SF%TMP_IGN) TW(IW) = T
   ENDIF
 
   ! Determine ghost wall temperature
 
   RHOWAL    = 0.5_EB*(RHO(IIG,JJG,KKG)+RHO_W(IW))
   CP_TERM   = MAX(0._EB,-CP_GAMMA*UW(IW)*RHOWAL)
   QCONF(IW) = HEAT_TRANS_COEF(IW) * (TMP_G - 0.5_EB * (TMP_F(IW) + TMP_F_OLD) )
   TMP_W(IW) = ( (RDN(IW)*KW(IW)-0.5_EB*CP_TERM)*TMP_G + CP_TERM*TMP_F(IW)-QCONF(IW) )/(0.5_EB*CP_TERM+RDN(IW)*KW(IW))
   TMP_W(IW) = MAX(TMPMIN,TMP_W(IW))
   IF (SOLID(CELL_INDEX(II,JJ,KK))) TMP(II,JJ,KK) = MAX(100._EB,MIN(4900._EB,TMP_W(IW)))

ENDDO WALL_CELL_LOOP
 
END SUBROUTINE PYROLYSIS
 
 
REAL(EB) FUNCTION HEAT_TRANSFER_COEFFICIENT(IW,IIG,JJG,KKG,IOR,TMP_G,DELTA_TMP)

INTEGER  :: IW,IIG,JJG,KKG,IOR,ITMP
REAL(EB) :: TMP_G,DELTA_TMP,U2,V2,W2,VELCON,H_NATURAL,H_FORCED,H_DNS
REAL(EB), POINTER, DIMENSION(:,:,:) :: UU,VV,WW,RHOP
 
IF (H_FIXED >= 0._EB) THEN
   HEAT_TRANSFER_COEFFICIENT = H_FIXED
   RETURN
ENDIF

IF (PREDICTOR) THEN
   UU => U
   VV => V
   WW => W
   RHOP => RHOS
ELSE
   UU => US
   VV => VS
   WW => WS
   RHOP => RHO
ENDIF 

IF (DNS) THEN
   HEAT_TRANSFER_COEFFICIENT = 2._EB*KW(IW)*RDN(IW)
ELSE
   IF (IIG<0) THEN   ! No gas cell information available
      SELECT CASE(ABS(IOR))
         CASE(1:2)
            H_NATURAL = HCV*ABS(DELTA_TMP)**ONTH
         CASE(3)
            H_NATURAL = HCH*ABS(DELTA_TMP)**ONTH
      END SELECT
      H_FORCED = 0._EB
   ELSE
      SELECT CASE(ABS(IOR))
         CASE(1)
            V2 = 0.25_EB*(VV(IIG,JJG,KKG)+VV(IIG,JJG-1,KKG))**2
            W2 = 0.25_EB*(WW(IIG,JJG,KKG)+WW(IIG,JJG,KKG-1))**2
            VELCON = (V2+W2)**0.4_EB
            H_NATURAL = HCV*ABS(DELTA_TMP)**ONTH
         CASE(2)
            U2 = 0.25_EB*(UU(IIG,JJG,KKG)+UU(IIG-1,JJG,KKG))**2
            W2 = 0.25_EB*(WW(IIG,JJG,KKG)+WW(IIG,JJG,KKG-1))**2
            VELCON = (U2+W2)**0.4_EB
            H_NATURAL = HCV*ABS(DELTA_TMP)**ONTH
         CASE(3)
            U2 = 0.25_EB*(UU(IIG,JJG,KKG)+UU(IIG-1,JJG,KKG))**2
            V2 = 0.25_EB*(VV(IIG,JJG,KKG)+VV(IIG,JJG-1,KKG))**2
            VELCON = (U2+V2)**0.4_EB
            H_NATURAL = HCH*ABS(DELTA_TMP)**ONTH
      END SELECT
      H_FORCED  = C_FORCED*VELCON*RHOP(IIG,JJG,KKG)**0.8_EB
   ENDIF
   ITMP = 0.1_EB*TMP_G
   H_DNS = SPECIES(0)%MU(ITMP)*CP_GAMMA*RPR*2._EB*RDN(IW)
   HEAT_TRANSFER_COEFFICIENT = MAX(H_DNS,H_FORCED,H_NATURAL)
ENDIF

END FUNCTION HEAT_TRANSFER_COEFFICIENT

END SUBROUTINE WALL_BC

SUBROUTINE GET_REV_wall(MODULE_REV,MODULE_DATE)
INTEGER,INTENT(INOUT) :: MODULE_REV
CHARACTER(255),INTENT(INOUT) :: MODULE_DATE

WRITE(MODULE_DATE,'(A)') wallrev(INDEX(wallrev,':')+1:LEN_TRIM(wallrev)-2)
READ (MODULE_DATE,'(I5)') MODULE_REV
WRITE(MODULE_DATE,'(A)') walldate

END SUBROUTINE GET_REV_wall
 
END MODULE WALL_ROUTINES