Merge pull request #14941 from cxp484/FireX

FireX: Merge with firemodels/master

Author: cxp484

Manuals/FDS_User_Guide/FDS_User_Guide.tex

@@ -464,6 +464,7 @@ MODULE GLOBAL_CONSTANTS
 INTEGER :: N_PASSIVE_SCALARS=0                                      !< Number of passive scalars
 INTEGER :: N_TOTAL_SCALARS=0                                        !< Number of total scalars, tracked and passive
 INTEGER :: N_FIXED_CHEMISTRY_SUBSTEPS=-1                            !< Number of chemistry substeps in combustion routine
+INTEGER :: ZETA_0_RAMP_INDEX=0                                      !< Ramp index for initial unmixed fraction
 
 LOGICAL :: OUTPUT_CHEM_IT=.FALSE.
 LOGICAL :: REAC_SOURCE_CHECK=.FALSE.

Source/cons.f90

@@ -96,12 +96,13 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT)
 USE PHYSICAL_FUNCTIONS, ONLY: GET_SPECIFIC_GAS_CONSTANT,GET_MASS_FRACTION_ALL,GET_SPECIFIC_HEAT,GET_MOLECULAR_WEIGHT, &
                               GET_SENSIBLE_ENTHALPY_Z,IS_REALIZABLE
 USE CHEMCONS, ONLY: DO_CHEM_LOAD_BALANCE
+USE MATH_FUNCTIONS, ONLY: EVALUATE_RAMP
 
 INTEGER :: NM,I,J,K,NC, CHEM_SUBIT_TMP,CHEM_SUBIT_TMP_OUT,ICC,JCC,INDX,NRECEIVE_CELLS
 INTEGER :: NCHEM_ACTIVE_CELLS_AND_CC, NCHEM_ACTIVE_CELLS_AND_CC_GLOBAL, IERR, NCHEM_ACTIVE_CELLS, NCHEM_ACTIVE_CC
 REAL(EB), INTENT(IN) :: T,DT
 REAL(EB) :: ZZ_GET(1:N_TRACKED_SPECIES),DZZ(1:N_TRACKED_SPECIES),&
-            REAC_SOURCE_TERM_TMP(N_TRACKED_SPECIES),Q_REAC_TMP(N_REACTIONS),VCELL,PRES
+            REAC_SOURCE_TERM_TMP(N_TRACKED_SPECIES),Q_REAC_TMP(N_REACTIONS),VCELL,PRES,ZETA_0,TSI
 LOGICAL :: Q_EXISTS
 TYPE (CC_CUTCELL_TYPE), POINTER :: CC
 LOGICAL :: IS_CHEM_ACTIVE
@@ -175,7 +176,10 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT)
 !       For 1 MPI process (serial) distribution is not needed.
 !------
 
- DO_CHEM_LOAD_BALANCE_IF : IF(N_MPI_PROCESSES > 1 .AND. DO_CHEM_LOAD_BALANCE) THEN
+TSI = T-T_BEGIN
+ZETA_0 = EVALUATE_RAMP(TSI,ZETA_0_RAMP_INDEX)*INITIAL_UNMIXED_FRACTION
+
+CHEM_LOAD_BALANCE_IF : IF(N_MPI_PROCESSES > 1 .AND. DO_CHEM_LOAD_BALANCE) THEN
    IF (.NOT. ALLOCATED(NCHEM_ACTIVE_CELLS_AND_CC_RANK_WISE)) ALLOCATE(NCHEM_ACTIVE_CELLS_AND_CC_RANK_WISE(N_MPI_PROCESSES))
    NCHEM_ACTIVE_CELLS_AND_CC_RANK_WISE = 0
    CALL MPI_ALLGATHER(NCHEM_ACTIVE_CELLS_AND_CC,1,MPI_INTEGER,NCHEM_ACTIVE_CELLS_AND_CC_RANK_WISE,1,MPI_INTEGER,MPI_COMM_WORLD,IERR)
@@ -205,7 +209,7 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT)
          ! Call combustion integration routine for Cartesian cell (I,J,K)
          CALL COMBUSTION_MODEL( T,DT,ZZ_GET,Q_OUT,MIX_TIME_OUT,CHI_R_OUT,&
                                 CHEM_SUBIT_TMP_OUT,REAC_SOURCE_TERM_TMP,Q_REAC_TMP,&
-                                MYTEMP,MYRHO,PRES,MYMU,DELTA,VOL)
+                                MYTEMP,MYRHO,PRES,MYMU,DELTA,VOL,ZETA_0)
          !***************************************************************************************
 
          RESULTS_TO_SEND_ARRAY(1:N_TRACKED_SPECIES,INDX) = ZZ_GET
@@ -216,7 +220,6 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT)
 
       ENDDO
 
-
       IF (COMBUSTION_ODE_SOLVER == CVODE_SOLVER) THEN
          T_CHEM_ODE = T_CHEM_ODE+CURRENT_TIME()-TNOW2
       ENDIF
@@ -230,7 +233,7 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT)
       IF (STOP_STATUS/=NO_STOP) RETURN
    ENDIF !NCHEM_ACTIVE_CELLS_AND_CC_GLOBAL >0
 
-ELSE DO_CHEM_LOAD_BALANCE_IF
+ELSE CHEM_LOAD_BALANCE_IF
 
    IF (NCHEM_ACTIVE_CELLS_AND_CC >0) THEN
       ! Serial chemistry:
@@ -254,7 +257,7 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT)
             CALL COMBUSTION_MODEL( T,DT,ZZ_GET,Q(I,J,K),MIX_TIME(I,J,K),CHI_R(I,J,K),&
                                    CHEM_SUBIT_TMP,REAC_SOURCE_TERM_TMP,Q_REAC_TMP,&
                                    TMP(I,J,K),RHO(I,J,K),PRES, MU(I,J,K),&
-                                   LES_FILTER_WIDTH(I,J,K),DX(I)*DY(J)*DZ(K),IIC=I,JJC=J,KKC=K )
+                                   LES_FILTER_WIDTH(I,J,K),DX(I)*DY(J)*DZ(K),ZETA_0,IIC=I,JJC=J,KKC=K )
             !***************************************************************************************
             !IF (STOP_STATUS/=NO_STOP) RETURN
             IF (OUTPUT_CHEM_IT) CHEM_SUBIT(I,J,K) = CHEM_SUBIT_TMP
@@ -278,7 +281,7 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT)
                                       CC%CHI_R(JCC),&
                                       CHEM_SUBIT_TMP,REAC_SOURCE_TERM_TMP,Q_REAC_TMP,&
                                       CC%TMP(JCC),CC%RHO(JCC),PRES,MU(I,J,K),&
-                                      LES_FILTER_WIDTH(I,J,K),CC%VOLUME(JCC),IIC=I,JJC=J,KKC=K)
+                                      LES_FILTER_WIDTH(I,J,K),CC%VOLUME(JCC),ZETA_0,IIC=I,JJC=J,KKC=K)
                !***************************************************************************************
                CALL SET_SPECIES_SOURCE_TERM_CUTCELL(DT, ICC, JCC, ZZ_GET, DZZ, REAC_SOURCE_TERM_TMP, Q_REAC_TMP)
             END DO
@@ -291,7 +294,7 @@ SUBROUTINE COMBUSTION_GENERAL_LOAD_BALANCED(T,DT)
          T_CHEM_ODE = T_CHEM_ODE+CURRENT_TIME()-TNOW2
       ENDIF
    ENDIF !NCHEM_ACTIVE_CELLS_AND_CC >0
-ENDIF DO_CHEM_LOAD_BALANCE_IF
+ENDIF CHEM_LOAD_BALANCE_IF
 
 ! This volume refactoring is needed for RADIATION_FVM (CHI_R, Q) and plotting slices:
 IF(CC_IBM) THEN
@@ -715,13 +718,13 @@ END SUBROUTINE GATHER_CELLS_FROM_MPI_PROCESSES
 
 
 SUBROUTINE COMBUSTION_MODEL(T,DT,ZZ_GET,Q_OUT,MIX_TIME_OUT,CHI_R_OUT,CHEM_SUBIT_OUT,REAC_SOURCE_TERM_OUT,Q_REAC_OUT,&
-                            TMP_IN,RHO_IN,PRES_IN,MU_IN,DELTA,CELL_VOLUME,IIC,JJC,KKC)
+                            TMP_IN,RHO_IN,PRES_IN,MU_IN,DELTA,CELL_VOLUME,ZETA_0_IN,IIC,JJC,KKC)
 USE MATH_FUNCTIONS, ONLY: EVALUATE_RAMP
 USE PHYSICAL_FUNCTIONS, ONLY: GET_REALIZABLE_MF
 USE COMP_FUNCTIONS, ONLY: SHUTDOWN
 USE CHEMCONS, ONLY: ODE_MIN_ATOL, FLAME_THICK_FACTOR
 INTEGER, INTENT(IN), OPTIONAL :: IIC,JJC,KKC
-REAL(EB), INTENT(IN) :: T,DT,RHO_IN,PRES_IN,MU_IN,DELTA,CELL_VOLUME
+REAL(EB), INTENT(IN) :: T,DT,RHO_IN,PRES_IN,MU_IN,DELTA,CELL_VOLUME,ZETA_0_IN
 REAL(EB), INTENT(OUT) :: Q_OUT,MIX_TIME_OUT,CHI_R_OUT,REAC_SOURCE_TERM_OUT(N_TRACKED_SPECIES),Q_REAC_OUT(N_REACTIONS)
 INTEGER, INTENT(OUT) :: CHEM_SUBIT_OUT
 REAL(EB), INTENT(INOUT) :: ZZ_GET(1:N_TRACKED_SPECIES)
@@ -733,7 +736,7 @@ SUBROUTINE COMBUSTION_MODEL(T,DT,ZZ_GET,Q_OUT,MIX_TIME_OUT,CHI_R_OUT,CHEM_SUBIT_
             Q_REAC_SUM(1:N_REACTIONS),Q_SUM_CHI_R,CHI_R_SUM,TIME_RAMP_FACTOR,&
             TOTAL_MIXED_MASS_1,TOTAL_MIXED_MASS_2,TOTAL_MIXED_MASS_4,TOTAL_MIXED_MASS,&
             ZETA_1,ZETA_2,ZETA_4,D_F,TMP_IN,C_U,DT_SUB_OLD,ERR_EST(N_TRACKED_SPECIES),ERR_TOL(N_TRACKED_SPECIES),ERR_TINY,&
-            ZZ_TEMP(1:N_TRACKED_SPECIES), ATOL(1:N_TRACKED_SPECIES)
+            ZZ_TEMP(1:N_TRACKED_SPECIES),ATOL(1:N_TRACKED_SPECIES)
 INTEGER :: NR,NS,ITER,TVI,RICH_ITER,TIME_ITER,RICH_ITER_MAX
 INTEGER, PARAMETER :: TV_ITER_MIN=5
 LOGICAL :: TV_FLUCT(1:N_TRACKED_SPECIES),EXTINCT,NO_REACTIONS,NO_REAC_2,NO_REAC_4
@@ -764,7 +767,7 @@ SUBROUTINE COMBUSTION_MODEL(T,DT,ZZ_GET,Q_OUT,MIX_TIME_OUT,CHI_R_OUT,CHEM_SUBIT_
    END SELECT
 ENDIF
 
-ZETA_0 = INITIAL_UNMIXED_FRACTION
+ZETA_0 = ZETA_0_IN
 CELL_MASS = RHO_IN*CELL_VOLUME
 
 DT_SUB_MIN = DT/REAL(MAX_CHEMISTRY_SUBSTEPS,EB)

Source/fire.f90

@@ -1849,7 +1849,7 @@ REAL(EB) FUNCTION GET_PARTICLE_ENTHALPY(I_LPC,TMP_S)
 USE MATH_FUNCTIONS, ONLY: INTERPOLATE1D_UNIFORM
 REAL(EB), INTENT(IN) :: TMP_S
 REAL(EB) :: RHO,RHO_H,VOL,DTMP,H_S,THICKNESS
-INTEGER :: I,N,ITMP,I_GRAD
+INTEGER :: I,N,ITMP
 INTEGER, INTENT(IN) :: I_LPC
 TYPE(LAGRANGIAN_PARTICLE_CLASS_TYPE), POINTER :: LPC
 TYPE(SURFACE_TYPE), POINTER :: SF
@@ -1861,21 +1861,16 @@ REAL(EB) FUNCTION GET_PARTICLE_ENTHALPY(I_LPC,TMP_S)
    CALL INTERPOLATE1D_UNIFORM(LBOUND(SPECIES(LPC%Y_INDEX)%H_L,1),SPECIES(LPC%Y_INDEX)%H_L,TMP_S,GET_PARTICLE_ENTHALPY)
 ELSE
    SF=>SURFACE(LPC%SURF_INDEX)
-   SELECT CASE(SF%GEOMETRY)
-      CASE(SURF_CARTESIAN)                          ; I_GRAD = 1
-      CASE(SURF_CYLINDRICAL,SURF_INNER_CYLINDRICAL) ; I_GRAD = 2
-      CASE(SURF_SPHERICAL)                          ; I_GRAD = 3
-   END SELECT
    RHO_H = 0._EB
    RHO = 0._EB
    ITMP = MIN(I_MAX_TEMP-1,INT(TMP_S))
    DTMP = TMP_S-REAL(ITMP,EB)
    THICKNESS = SUM(SF%LAYER_THICKNESS)
    DO I=1,SUM(SF%N_LAYER_CELLS)
       IF (SF%GEOMETRY==SURF_INNER_CYLINDRICAL) THEN
-         VOL = (SF%INNER_RADIUS+SF%X_S(I))**I_GRAD - (SF%INNER_RADIUS+SF%X_S(I-1))**I_GRAD
+         VOL = (SF%INNER_RADIUS+SF%X_S(I))**SF%I_GRAD - (SF%INNER_RADIUS+SF%X_S(I-1))**SF%I_GRAD
       ELSE
-         VOL = (THICKNESS+SF%INNER_RADIUS-SF%X_S(I-1))**I_GRAD - (THICKNESS+SF%INNER_RADIUS-SF%X_S(I))**I_GRAD
+         VOL = (THICKNESS+SF%INNER_RADIUS-SF%X_S(I-1))**SF%I_GRAD - (THICKNESS+SF%INNER_RADIUS-SF%X_S(I))**SF%I_GRAD
       ENDIF
       MATL_REMESH: DO N=1,SF%N_MATL
          IF (SF%RHO_0(I,N)<=TWO_EPSILON_EB) CYCLE MATL_REMESH
@@ -2121,7 +2116,7 @@ REAL(EB) FUNCTION SURFACE_DENSITY(MODE,SF,ONE_D,MATL_INDEX)
 
 INTEGER, INTENT(IN) :: MODE
 INTEGER, INTENT(IN), OPTIONAL :: MATL_INDEX
-INTEGER :: I_GRAD,NWP,II2,N,ITMP
+INTEGER :: NWP,II2,N,ITMP
 REAL(EB) :: WGT,R_S(0:NWP_MAX),EPUM,DTMP
 TYPE(BOUNDARY_ONE_D_TYPE), INTENT(IN), POINTER :: ONE_D
 TYPE(SURFACE_TYPE), INTENT(IN), POINTER :: SF
@@ -2133,12 +2128,6 @@ REAL(EB) FUNCTION SURFACE_DENSITY(MODE,SF,ONE_D,MATL_INDEX)
 
 ELSE THERMALLY_THICK_IF
 
-   SELECT CASE(SF%GEOMETRY)
-      CASE(SURF_CARTESIAN)                           ; I_GRAD = 1
-      CASE(SURF_CYLINDRICAL,SURF_INNER_CYLINDRICAL)  ; I_GRAD = 2
-      CASE(SURF_SPHERICAL)                           ; I_GRAD = 3
-   END SELECT
-
    NWP = SUM(ONE_D%N_LAYER_CELLS)
    IF (SF%GEOMETRY==SURF_INNER_CYLINDRICAL) THEN
       R_S(0:NWP) = SF%INNER_RADIUS + ONE_D%X(0:NWP)
@@ -2149,8 +2138,9 @@ REAL(EB) FUNCTION SURFACE_DENSITY(MODE,SF,ONE_D,MATL_INDEX)
    SURFACE_DENSITY = 0._EB
    NUMBER_WALL_POINTS_LOOP: DO II2=1,NWP
       AREA_VOLUME_SELECT: SELECT CASE(MODE)
-         CASE(0,2); WGT = ABS(R_S(II2-1)**I_GRAD-R_S(II2)**I_GRAD)/(REAL(I_GRAD,EB)*(SF%INNER_RADIUS+SF%THICKNESS)**(I_GRAD-1))
-         CASE(1,3); WGT = ABS(R_S(II2-1)**I_GRAD-R_S(II2)**I_GRAD)/(SF%INNER_RADIUS+SF%THICKNESS)**I_GRAD
+         CASE(0,2); WGT = ABS(R_S(II2-1)**SF%I_GRAD-R_S(II2)**SF%I_GRAD)/ &
+                          (REAL(SF%I_GRAD,EB)*(SF%INNER_RADIUS+SF%THICKNESS)**(SF%I_GRAD-1))
+         CASE(1,3); WGT = ABS(R_S(II2-1)**SF%I_GRAD-R_S(II2)**SF%I_GRAD)/(SF%INNER_RADIUS+SF%THICKNESS)**SF%I_GRAD
       END SELECT AREA_VOLUME_SELECT
 
       EPUM = 1._EB ! energy per unit mass
@@ -4072,6 +4062,7 @@ SUBROUTINE PACK_BOUNDARY_ONE_D(NM,IC,RC,LC,OS,OD_INDEX,UNPACK_IT,COUNT_ONLY,CHEC
 I1=RC+1 ; RC=I1+ONE_D%N_LAYERS-1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(I1:RC),ONE_D%LAYER_THICKNESS(1:RC-I1+1)    , UNPACK_IT)
 I1=RC+1 ; RC=I1+ONE_D%N_LAYERS-1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(I1:RC),ONE_D%LAYER_THICKNESS_OLD(1:RC-I1+1), UNPACK_IT)
 I1=RC+1 ; RC=I1+ONE_D%N_LAYERS-1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(I1:RC),ONE_D%MIN_LAYER_THICKNESS(1:RC-I1+1), UNPACK_IT)
+I1=RC+1 ; RC=I1+ONE_D%N_LAYERS-1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(I1:RC),ONE_D%MIN_LAYER_MASS(1:RC-I1+1), UNPACK_IT)
 I1=RC+1 ; RC=I1+ONE_D%N_LAYERS-1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(I1:RC),ONE_D%MIN_DIFFUSIVITY(1:RC-I1+1)    , UNPACK_IT)
 I1=RC+1 ; RC=I1+ONE_D%N_LAYERS-1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(I1:RC),ONE_D%DDSUM(1:RC-I1+1)              , UNPACK_IT)
 I1=RC+1 ; RC=I1+ONE_D%N_LAYERS-1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(I1:RC),ONE_D%SMALLEST_CELL_SIZE(1:RC-I1+1) , UNPACK_IT)
@@ -4115,6 +4106,7 @@ SUBROUTINE REALLOCATE_BOUNDARY_ONE_D(ONE_D)
    ALLOCATE(ONE_D%LAYER_THICKNESS_OLD(ONE_D%N_LAYERS))
 IF (ALLOCATED(ONE_D%MIN_LAYER_THICKNESS)) DEALLOCATE(ONE_D%MIN_LAYER_THICKNESS)
    ALLOCATE(ONE_D%MIN_LAYER_THICKNESS(ONE_D%N_LAYERS))
+IF (ALLOCATED(ONE_D%MIN_LAYER_MASS)) DEALLOCATE(ONE_D%MIN_LAYER_MASS)         ; ALLOCATE(ONE_D%MIN_LAYER_MASS(ONE_D%N_LAYERS))
 IF (ALLOCATED(ONE_D%HT3D_LAYER)) DEALLOCATE(ONE_D%HT3D_LAYER)                 ; ALLOCATE(ONE_D%HT3D_LAYER(ONE_D%N_LAYERS))
 IF (ALLOCATED(ONE_D%MIN_DIFFUSIVITY)) DEALLOCATE(ONE_D%MIN_DIFFUSIVITY)       ; ALLOCATE(ONE_D%MIN_DIFFUSIVITY(ONE_D%N_LAYERS))
 IF (ALLOCATED(ONE_D%RHO_C_S)) DEALLOCATE(ONE_D%RHO_C_S)                       ; ALLOCATE(ONE_D%RHO_C_S(ONE_D%N_CELLS_MAX))
@@ -4141,12 +4133,12 @@ END SUBROUTINE REALLOCATE_BOUNDARY_ONE_D
 
 SUBROUTINE INITIALIZE_BOUNDARY_ONE_D(NM,OD_INDEX,SURF_INDEX)
 
-USE GLOBAL_CONSTANTS, ONLY: RADIATION
+USE GLOBAL_CONSTANTS, ONLY: RADIATION,NWP_MAX
 INTEGER, INTENT(IN) :: NM,OD_INDEX,SURF_INDEX
 TYPE(BOUNDARY_ONE_D_TYPE), POINTER :: ONE_D
 TYPE(SURFACE_TYPE), POINTER :: SF
 INTEGER :: NN,I,II
-REAL(EB) :: RAMP_POSITION,TF,TB
+REAL(EB) :: RAMP_POSITION,TF,TB,R(0:NWP_MAX)
 
 ONE_D => MESHES(NM)%BOUNDARY_ONE_D(OD_INDEX)
 SF => SURFACE(SURF_INDEX)
@@ -4162,6 +4154,12 @@ SUBROUTINE INITIALIZE_BOUNDARY_ONE_D(NM,OD_INDEX,SURF_INDEX)
 DO I=1,MIN(ONE_D%N_CELLS_MAX,ONE_D%N_CELLS_INI)
    ONE_D%DX_OLD(I) = ONE_D%X(I)-ONE_D%X(I-1)
 ENDDO
+R(0:ONE_D%N_CELLS_INI) = SF%THICKNESS + SF%INNER_RADIUS - ONE_D%X(0:ONE_D%N_CELLS_INI)
+I = 0
+DO NN = 1,ONE_D%N_LAYERS
+   ONE_D%MIN_LAYER_MASS(NN) = SF%MIN_LAYER_MASS(NN)*(R(I)**SF%I_GRAD - R(I+ONE_D%N_LAYER_CELLS(NN))**SF%I_GRAD)
+   I = I + ONE_D%N_LAYER_CELLS(NN)
+ENDDO
 IF (ALLOCATED(ONE_D%LAYER_THICKNESS_OLD)) ONE_D%LAYER_THICKNESS_OLD(1:ONE_D%N_LAYERS) = SF%LAYER_THICKNESS(1:SF%N_LAYERS)
 IF (SF%RAMP_T_I_INDEX > 0) THEN
    !NOTE: Replicating EVALUATE_RAMP since MODULE MATH_FUNCTIONS uses the MODULE which contains this routine
@@ -4319,6 +4317,7 @@ SUBROUTINE PACK_BOUNDARY_PROP1(NM,IC,RC,LC,OS,B1_INDEX,UNPACK_IT,COUNT_ONLY,SURF
 RC=RC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(RC),B1%RHO_G,UNPACK_IT)
 RC=RC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(RC),B1%RDN,UNPACK_IT)
 RC=RC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(RC),B1%K_G,UNPACK_IT)
+RC=RC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(RC),B1%M_DOT_LAYER_PP,UNPACK_IT)
 RC=RC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(RC),B1%Q_DOT_G_PP,UNPACK_IT)
 RC=RC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(RC),B1%Q_DOT_O2_PP,UNPACK_IT)
 RC=RC+1 ; IF (.NOT.COUNT_ONLY) CALL EQUATE(OS%REALS(RC),B1%Q_CONDENSE,UNPACK_IT)
@@ -5267,24 +5266,14 @@ END SUBROUTINE GET_WALL_NODE_WEIGHTS
 !> \param X_S_NEW Array of interior cell edge positions after shrinkage or swelling (m)
 !> \param INT_WGT Array of weighting factors for new arrangement of interior cells
 
-SUBROUTINE GET_INTERPOLATION_WEIGHTS(GEOMETRY,NWP,NWP_NEW,INNER_RADIUS,X_S,X_S_NEW,INT_WGT)
+SUBROUTINE GET_INTERPOLATION_WEIGHTS(GEOMETRY,I_GRAD,NWP,NWP_NEW,INNER_RADIUS,X_S,X_S_NEW,INT_WGT)
 
-INTEGER, INTENT(IN)  :: GEOMETRY,NWP,NWP_NEW
+INTEGER, INTENT(IN)  :: GEOMETRY,I_GRAD,NWP,NWP_NEW
 REAL(EB), INTENT(IN) :: X_S(0:NWP), X_S_NEW(0:NWP_NEW), INNER_RADIUS
 REAL(EB), INTENT(OUT) :: INT_WGT(NWP_NEW,NWP)
 
 REAL(EB) XUP,XLOW,XUP_NEW,XLOW_NEW,VOL_NEW,VOL,THICKNESS
-INTEGER I_NEW, I_OLD, I_GRAD
-
-
-SELECT CASE(GEOMETRY)
-   CASE(SURF_CARTESIAN)
-      I_GRAD = 1
-   CASE(SURF_CYLINDRICAL,SURF_INNER_CYLINDRICAL)
-      I_GRAD = 2
-   CASE(SURF_SPHERICAL)
-      I_GRAD = 3
-END SELECT
+INTEGER I_NEW, I_OLD
 
 I_OLD = 1
 INT_WGT = 0._EB

Source/func.f90

@@ -4188,6 +4188,7 @@ SUBROUTINE FIND_WALL_BACK_INDEX(NM,IW)
    DEALLOCATE(ONE_D%MATL_INDEX)      ; ALLOCATE(ONE_D%MATL_INDEX(ONE_D%N_MATL))
    DEALLOCATE(ONE_D%LAYER_THICKNESS) ; ALLOCATE(ONE_D%LAYER_THICKNESS(ONE_D%N_LAYERS))
    DEALLOCATE(ONE_D%MIN_LAYER_THICKNESS) ; ALLOCATE(ONE_D%MIN_LAYER_THICKNESS(ONE_D%N_LAYERS))
+   DEALLOCATE(ONE_D%MIN_LAYER_MASS) ; ALLOCATE(ONE_D%MIN_LAYER_MASS(ONE_D%N_LAYERS))
    DEALLOCATE(ONE_D%HT3D_LAYER)      ; ALLOCATE(ONE_D%HT3D_LAYER(ONE_D%N_LAYERS))
    ONE_D%LAYER_THICKNESS(1:ONE_D%N_LAYERS) = LAYER_THICKNESS(1:ONE_D%N_LAYERS)
    ONE_D%MIN_LAYER_THICKNESS(1:ONE_D%N_LAYERS) = MIN_LAYER_THICKNESS(1:ONE_D%N_LAYERS)