Merge pull request #13889 from rmcdermo/master

FDS Source: update to FLUX_LIMITER_MW_CORRECTION

Author: rmcdermo

Source/ccib.f90

@@ -8130,13 +8130,13 @@ SUBROUTINE SET_EXIMRHOZZLIM_3D(NM,N)
 
 ! Local Variables:
 INTEGER :: I,J,K,X1AXIS
-REAL(EB), POINTER, DIMENSION(:,:,:) :: FX_ZZ=>NULL(),FY_ZZ=>NULL(),FZ_ZZ=>NULL()
+REAL(EB), POINTER, DIMENSION(:,:,:,:) :: FX_ZZ=>NULL(),FY_ZZ=>NULL(),FZ_ZZ=>NULL()
 INTEGER :: IFACE,IW
 TYPE(CC_REGFACEZ_TYPE),  POINTER, DIMENSION(:) :: REGFACE_Z=>NULL()
 
-FX_ZZ=>WORK2
-FY_ZZ=>WORK3
-FZ_ZZ=>WORK4
+FX_ZZ=>SWORK1
+FY_ZZ=>SWORK2
+FZ_ZZ=>SWORK3
 
 AXIS_DO : DO X1AXIS = IAXIS,KAXIS
    SELECT CASE(X1AXIS)
@@ -8164,11 +8164,11 @@ SUBROUTINE SET_EXIMRHOZZLIM_3D(NM,N)
       ELSE
          SELECT CASE(X1AXIS)
          CASE(IAXIS)
-            REGFACE_Z(IFACE)%FN_ZZ(N) = FX_ZZ(I,J,K)
+            REGFACE_Z(IFACE)%FN_ZZ(N) = FX_ZZ(I,J,K,N)
          CASE(JAXIS)
-            REGFACE_Z(IFACE)%FN_ZZ(N) = FY_ZZ(I,J,K)
+            REGFACE_Z(IFACE)%FN_ZZ(N) = FY_ZZ(I,J,K,N)
          CASE(KAXIS)
-            REGFACE_Z(IFACE)%FN_ZZ(N) = FZ_ZZ(I,J,K)
+            REGFACE_Z(IFACE)%FN_ZZ(N) = FZ_ZZ(I,J,K,N)
          END SELECT
       ENDIF
    ENDDO IFACE_DO

Source/cons.f90

@@ -270,7 +270,7 @@ MODULE GLOBAL_CONSTANTS
 LOGICAL :: TENSOR_DIFFUSIVITY=.FALSE.               !< If true, use experimental tensor diffusivity model for spec and tmp
 LOGICAL :: OXPYRO_MODEL=.FALSE.                     !< Flag to use oxidative pyrolysis mass transfer model
 LOGICAL :: OUTPUT_WALL_QUANTITIES=.FALSE.           !< Flag to force call to WALL_MODEL
-LOGICAL :: TEST_FLUX_LIMITER_FACE_CORRECTION=.FALSE.
+LOGICAL :: FLUX_LIMITER_MW_CORRECTION=.FALSE.
 
 INTEGER, ALLOCATABLE, DIMENSION(:) :: CHANGE_TIME_STEP_INDEX      !< Flag to indicate if a mesh needs to change time step
 INTEGER, ALLOCATABLE, DIMENSION(:) :: SETUP_PRESSURE_ZONES_INDEX  !< Flag to indicate if a mesh needs to keep searching for ZONEs

Source/divg.f90

@@ -1616,7 +1616,7 @@ SUBROUTINE CONDENSATION_EVAPORATION(DT,NM)
 
 CALL POINT_TO_MESH(NM)
 
-ZZ_INTERIM=> SCALAR_WORK1
+ZZ_INTERIM=> SWORK1
 ZZ_INTERIM = ZZ
 RHO_INTERIM => WORK1
 RHO_INTERIM = RHO

Source/fire.f90

@@ -188,18 +188,15 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM)
 ! Allocate scalar face values
 
 ! Required for cell face density correction for multicomponent mixtures
-IF (TEST_FLUX_LIMITER_FACE_CORRECTION .AND. N_TRACKED_SPECIES>2) THEN
+IF (FLUX_LIMITER_MW_CORRECTION) THEN
    N_LOWER_SCALARS=0
 ELSE
    N_LOWER_SCALARS=1
 ENDIF
 
-ALLOCATE( M%FX(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO)
-CALL ChkMemErr('INIT','FX',IZERO)
-ALLOCATE( M%FY(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO)
-CALL ChkMemErr('INIT','FY',IZERO)
-ALLOCATE( M%FZ(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO)
-CALL ChkMemErr('INIT','FZ',IZERO)
+ALLOCATE( M%FX(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FX',IZERO)
+ALLOCATE( M%FY(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FY',IZERO)
+ALLOCATE( M%FZ(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FZ',IZERO)
 M%FX = 0._EB
 M%FY = 0._EB
 M%FZ = 0._EB
@@ -306,16 +303,16 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM)
 ALLOCATE(M%WORK8(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO) ; CALL ChkMemErr('INIT','WORK8',IZERO)
 ALLOCATE(M%WORK9(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO) ; CALL ChkMemErr('INIT','WORK9',IZERO)
 
-ALLOCATE(M%IWORK1(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO)                         ; CALL ChkMemErr('INIT','IWORK1',IZERO)
-ALLOCATE(M%SCALAR_WORK1(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO) ; CALL ChkMemErr('INIT','SCALAR_WORK1',IZERO)
-ALLOCATE(M%SCALAR_WORK2(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO) ; CALL ChkMemErr('INIT','SCALAR_WORK2',IZERO)
-ALLOCATE(M%SCALAR_WORK3(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO) ; CALL ChkMemErr('INIT','SCALAR_WORK3',IZERO)
-ALLOCATE(M%SCALAR_WORK4(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO) ; CALL ChkMemErr('INIT','SCALAR_WORK4',IZERO)
+ALLOCATE(M%IWORK1(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO)                                 ; CALL ChkMemErr('INIT','IWORK1',IZERO)
+ALLOCATE(M%SWORK1(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO) ; CALL ChkMemErr('INIT','SWORK1',IZERO)
+ALLOCATE(M%SWORK2(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO) ; CALL ChkMemErr('INIT','SWORK2',IZERO)
+ALLOCATE(M%SWORK3(0:IBP1,0:JBP1,0:KBP1,N_LOWER_SCALARS:N_TOTAL_SCALARS),STAT=IZERO) ; CALL ChkMemErr('INIT','SWORK3',IZERO)
+ALLOCATE(M%SWORK4(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO)               ; CALL ChkMemErr('INIT','SWORK4',IZERO)
 M%IWORK1=0
-M%SCALAR_WORK1=0._EB
-M%SCALAR_WORK2=0._EB
-M%SCALAR_WORK3=0._EB
-M%SCALAR_WORK4=0._EB
+M%SWORK1=0._EB
+M%SWORK2=0._EB
+M%SWORK3=0._EB
+M%SWORK4=0._EB
 
 IF (STRATIFICATION) THEN
 

Source/init.f90

@@ -178,7 +178,7 @@ SUBROUTINE MASS_FINITE_DIFFERENCES(NM)
 
 ENDDO SPECIES_LOOP
 
-FACE_CORRECTION_IF: IF (TEST_FLUX_LIMITER_FACE_CORRECTION) THEN
+FACE_CORRECTION_IF: IF (FLUX_LIMITER_MW_CORRECTION) THEN
 
    ! Repeat the above for DENSITY
 
@@ -375,7 +375,7 @@ SUBROUTINE DENSITY(T,DT,NM)
 UU=>WORK1
 VV=>WORK2
 WW=>WORK3
-DEL_RHO_D_DEL_Z__0=>SCALAR_WORK4
+DEL_RHO_D_DEL_Z__0=>SWORK4
 
 PREDICTOR_STEP: SELECT CASE (PREDICTOR)
 

Source/mass.f90

@@ -107,10 +107,10 @@ MODULE MESH_VARIABLES
 
    ! Work arrays
 
-   REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: SCALAR_WORK1
-   REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: SCALAR_WORK2
-   REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: SCALAR_WORK3
-   REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: SCALAR_WORK4
+   REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: SWORK1
+   REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: SWORK2
+   REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: SWORK3
+   REAL(EB), ALLOCATABLE, DIMENSION(:,:,:,:) :: SWORK4
    REAL(EB), ALLOCATABLE, DIMENSION(:,:,:) :: WORK1,WORK2,WORK3,WORK4,WORK5,WORK6,WORK7,WORK8,WORK9
    INTEGER, ALLOCATABLE, DIMENSION(:,:,:) :: IWORK1
    REAL(EB),     ALLOCATABLE, DIMENSION(:,:,:) :: PWORK1,PWORK2,PWORK3,PWORK4
@@ -349,7 +349,7 @@ MODULE MESH_POINTERS
    MU,MU_DNS,TMP,Q,KAPPA_GAS,CHI_R,QR,QR_W,UII,RSUM,D_SOURCE, &
    CSD2,MTR,MSR,WEM,MIX_TIME,CHEM_SUBIT,STRAIN_RATE,D_Z_MAX,PP_RESIDUAL,LES_FILTER_WIDTH
 REAL(EB), POINTER, DIMENSION(:,:,:,:) :: ZZ,ZZS,REAC_SOURCE_TERM,DEL_RHO_D_DEL_Z,FX,FY,FZ, &
-                                         SCALAR_WORK1,SCALAR_WORK2,SCALAR_WORK3,SCALAR_WORK4, &
+                                         SWORK1,SWORK2,SWORK3,SWORK4, &
                                          Q_REAC,AVG_DROP_DEN,AVG_DROP_TMP,AVG_DROP_RAD,AVG_DROP_AREA,M_DOT_PPP, &
                                          ADV_FX,ADV_FY,ADV_FZ,DIF_FX,DIF_FY,DIF_FZ,DIF_FXS,DIF_FYS,DIF_FZS
 INTEGER, POINTER, DIMENSION(:,:)   :: CHEM_ACTIVE_CELLS
@@ -551,10 +551,10 @@ SUBROUTINE POINT_TO_MESH(NM)
 DIF_FXS=>M%DIF_FXS
 DIF_FYS=>M%DIF_FYS
 DIF_FZS=>M%DIF_FZS
-SCALAR_WORK1=>M%SCALAR_WORK1
-SCALAR_WORK2=>M%SCALAR_WORK2
-SCALAR_WORK3=>M%SCALAR_WORK3
-SCALAR_WORK4=>M%SCALAR_WORK4
+SWORK1=>M%SWORK1
+SWORK2=>M%SWORK2
+SWORK3=>M%SWORK3
+SWORK4=>M%SWORK4
 WORK=>M%WORK
 LSAVE=>M%LSAVE
 LWORK=>M%LWORK

Source/mesh.f90

@@ -3634,7 +3634,7 @@ SUBROUTINE PARTICLE_MASS_ENERGY_TRANSFER(T,DT,NM)
 
    RHO_INTERIM => WORK1 ; RHO_INTERIM = RHO
    TMP_INTERIM => WORK2 ; TMP_INTERIM = TMP
-   ZZ_INTERIM  => SCALAR_WORK1 ; ZZ_INTERIM = ZZ
+   ZZ_INTERIM  => SWORK1 ; ZZ_INTERIM = ZZ
 
 ENDIF
 

Source/part.f90

@@ -1745,7 +1745,7 @@ SUBROUTINE READ_MISC
                 RAMP_UX,RAMP_UY,RAMP_UZ,RAMP_VX,RAMP_VY,RAMP_VZ,RAMP_WX,RAMP_WY,RAMP_WZ,&
                 RESTART,RESTART_CHID,SC,&
                 RND_SEED,SIMULATION_MODE,SMOKE3D_16,SMOKE_ALBEDO,SOLID_PHASE_ONLY,SOOT_DENSITY,SOOT_OXIDATION,&
-                TAU_DEFAULT,TENSOR_DIFFUSIVITY,TERRAIN_IMAGE,TEST_FLUX_LIMITER_FACE_CORRECTION,TEXTURE_ORIGIN,&
+                TAU_DEFAULT,TENSOR_DIFFUSIVITY,TERRAIN_IMAGE,FLUX_LIMITER_MW_CORRECTION,TEXTURE_ORIGIN,&
                 THERMOPHORETIC_DEPOSITION,THERMOPHORETIC_SETTLING,THICKEN_OBSTRUCTIONS,&
                 TMPA,TURBULENCE_MODEL,TURBULENT_DEPOSITION,UVW_FILE,&
                 VERBOSE,VISIBILITY_FACTOR,VN_MAX,VN_MIN,Y_CO2_INFTY,Y_O2_INFTY,&
@@ -2823,7 +2823,12 @@ SUBROUTINE READ_SPEC
 N_TRACKED_SPECIES = N_LUMPED + N_COPY + N_TOTAL_BINS + N_PREDEFINED_SMIX + N_CONDENSATION
 N_TOTAL_SCALARS = N_TRACKED_SPECIES + N_PASSIVE_SCALARS
 
-! Allocate the primitive species array.
+! Only allow use of flux limiter molecular weight correction if N_TRACKED_SPECIES>2
+
+IF (N_TRACKED_SPECIES<=2) FLUX_LIMITER_MW_CORRECTION=.FALSE.
+
+! Allocate the primitive species array
+
 ALLOCATE(SPECIES(N_SPECIES),STAT=IZERO)
 CALL ChkMemErr('READ','SPECIES',IZERO)
 ALLOCATE(Y_INDEX(N_SPECIES))

Source/read.f90

@@ -1651,9 +1651,9 @@ SUBROUTINE TENSOR_DIFFUSIVITY_MODEL(NM,OPT_N)
    ! SGS scalar flux
    ! CAUTION: The flux arrays must point to the same work arrays used in DIVERGENCE_PART_1
    ! Note: Do not reinitialize!  RHO_D_DZDX, etc., already store molecular diffusive flux
-   RHO_D_DZDX=>SCALAR_WORK1
-   RHO_D_DZDY=>SCALAR_WORK2
-   RHO_D_DZDZ=>SCALAR_WORK3
+   RHO_D_DZDX=>SWORK1
+   RHO_D_DZDY=>SWORK2
+   RHO_D_DZDZ=>SWORK3
 
    IF (PREDICTOR) THEN
       UU=>U