Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 0 additions & 1 deletion Source/cons.f90
Original file line number Diff line number Diff line change
Expand Up @@ -274,7 +274,6 @@ 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 :: FLUX_LIMITER_SINGLE_COEF=.TRUE. !< Flag to base flux coefficients off of a single worst-case scalar gradient
LOGICAL :: STORE_FIRE_ARRIVAL=.FALSE. !< Flag for tracking arrival of spreading fire front over a surface
LOGICAL :: STORE_FIRE_RESIDENCE=.FALSE. !< Flag for tracking residence time of spreading fire front over a surface
LOGICAL :: STORE_LS_SPREAD_RATE=.FALSE. !< Flag for outputting local level set spread rate magnitude
Expand Down
158 changes: 29 additions & 129 deletions Source/divg.f90
Original file line number Diff line number Diff line change
Expand Up @@ -586,7 +586,7 @@ SUBROUTINE DIVERGENCE_PART_1(T,DT,NM)

CONST_GAMMA_IF_1: IF (.NOT.CONSTANT_SPECIFIC_HEAT_RATIO) THEN

CALL ENTHALPY_ADVECTION(U_DOT_DEL_RHO_H_S) ! Compute u dot grad rho h_s
CALL ENTHALPY_ADVECTION_NEW(U_DOT_DEL_RHO_H_S) ! Compute u dot grad rho h_s

DO K=1,KBAR
DO J=1,JBAR
Expand Down Expand Up @@ -635,11 +635,7 @@ SUBROUTINE DIVERGENCE_PART_1(T,DT,NM)

CONST_GAMMA_IF_2: IF (.NOT.CONSTANT_SPECIFIC_HEAT_RATIO) THEN

IF (FLUX_LIMITER_SINGLE_COEF) THEN
CALL SPECIES_ADVECTION_PART_1_NEW
ELSE
CALL SPECIES_ADVECTION_PART_1 ! Compute and store face values of (rho Z_n)
ENDIF
CALL SPECIES_ADVECTION_PART_1_NEW ! Compute and store face values of (rho Z_n)

DO N=1,N_TRACKED_SPECIES

Expand Down Expand Up @@ -789,14 +785,14 @@ SUBROUTINE DIVERGENCE_PART_1(T,DT,NM)
END SUBROUTINE DIVERGENCE_PART_1


SUBROUTINE ENTHALPY_ADVECTION(U_DOT_DEL_RHO_H_S)
SUBROUTINE ENTHALPY_ADVECTION_NEW(U_DOT_DEL_RHO_H_S)

USE PHYSICAL_FUNCTIONS, ONLY: GET_SENSIBLE_ENTHALPY
USE MATH_FUNCTIONS, ONLY: GET_SCALAR_FACE_VALUE
USE MATH_FUNCTIONS, ONLY: GET_SCALAR_FACE_COEF,GET_SCALAR_FACE_VALUE_NEW
REAL(EB), POINTER, DIMENSION(:,:,:) :: FX_H_S,FY_H_S,FZ_H_S,RHO_H_S_P,U_DOT_DEL_RHO_H_S
REAL(EB) :: UN,UN_P,TMP_F_GAS,DU_P,DU_M,DV_P,DV_M,DW_P,DW_M,DU,H_S
REAL(EB), ALLOCATABLE, DIMENSION(:) :: ZZ_GET
REAL(EB), DIMENSION(0:3,0:3,0:3) :: U_TEMP,F_TEMP
REAL(EB), DIMENSION(0:3,0:3,0:3) :: U_TEMP,F_TEMP,B_TEMP
REAL(EB), DIMENSION(-1:3,-1:3,-1:3) :: Z_TEMP
INTEGER :: IC,I,J,K,IW
TYPE(WALL_TYPE), POINTER :: WC
Expand All @@ -809,6 +805,10 @@ SUBROUTINE ENTHALPY_ADVECTION(U_DOT_DEL_RHO_H_S)
FZ_H_S=>WORK4 ; FZ_H_S = 0._EB
U_DOT_DEL_RHO_H_S=>WORK6 ; U_DOT_DEL_RHO_H_S=0._EB

BFX = 1.E10_EB
BFY = 1.E10_EB
BFZ = 1.E10_EB

! Compute and store rho*h_s

!$OMP PARALLEL PRIVATE(ZZ_GET, H_S)
Expand All @@ -829,9 +829,13 @@ SUBROUTINE ENTHALPY_ADVECTION(U_DOT_DEL_RHO_H_S)

! Compute scalar face values

CALL GET_SCALAR_FACE_VALUE(UU,RHO_H_S_P,FX_H_S,1,IBM1,1,JBAR,1,KBAR,1,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_VALUE(VV,RHO_H_S_P,FY_H_S,1,IBAR,1,JBM1,1,KBAR,2,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_VALUE(WW,RHO_H_S_P,FZ_H_S,1,IBAR,1,JBAR,1,KBM1,3,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_COEF(UU,RHO_H_S_P,BFX,0,IBAR,1,JBAR,1,KBAR,1,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_COEF(VV,RHO_H_S_P,BFY,1,IBAR,0,JBAR,1,KBAR,2,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_COEF(WW,RHO_H_S_P,BFZ,1,IBAR,1,JBAR,0,KBAR,3,I_FLUX_LIMITER)

CALL GET_SCALAR_FACE_VALUE_NEW(UU,RHO_H_S_P,FX_H_S,BFX,0,IBAR,1,JBAR,1,KBAR,1,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_VALUE_NEW(VV,RHO_H_S_P,FY_H_S,BFY,1,IBAR,0,JBAR,1,KBAR,2,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_VALUE_NEW(WW,RHO_H_S_P,FZ_H_S,BFZ,1,IBAR,1,JBAR,0,KBAR,3,I_FLUX_LIMITER)

ALLOCATE(ZZ_GET(1:N_TRACKED_SPECIES))

Expand Down Expand Up @@ -871,7 +875,8 @@ SUBROUTINE ENTHALPY_ADVECTION(U_DOT_DEL_RHO_H_S)
IF ((UU(BC%II+1,BC%JJ,BC%KK)>0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II+1,BC%JJ,BC%KK))%WALL_INDEX(+1)>0)) THEN
Z_TEMP(0:2,1,1) = (/RHO_H_S_P(BC%II+1,BC%JJ,BC%KK),RHO_H_S_P(BC%II+1:BC%II+2,BC%JJ,BC%KK)/)
U_TEMP(1,1,1) = UU(BC%II+1,BC%JJ,BC%KK)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,1,I_FLUX_LIMITER)
B_TEMP(1,1,1) = BFX(BC%II+1,BC%JJ,BC%KK)
CALL GET_SCALAR_FACE_VALUE_NEW(U_TEMP,Z_TEMP,F_TEMP,B_TEMP,1,1,1,1,1,1,1,I_FLUX_LIMITER)
FX_H_S(BC%II+1,BC%JJ,BC%KK) = F_TEMP(1,1,1)
ENDIF
CASE(-1) OFF_WALL_SELECT_1
Expand All @@ -881,35 +886,40 @@ SUBROUTINE ENTHALPY_ADVECTION(U_DOT_DEL_RHO_H_S)
IF ((UU(BC%II-2,BC%JJ,BC%KK)<0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II-1,BC%JJ,BC%KK))%WALL_INDEX(-1)>0)) THEN
Z_TEMP(1:3,1,1) = (/RHO_H_S_P(BC%II-2:BC%II-1,BC%JJ,BC%KK),RHO_H_S_P(BC%II-1,BC%JJ,BC%KK)/)
U_TEMP(1,1,1) = UU(BC%II-2,BC%JJ,BC%KK)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,1,I_FLUX_LIMITER)
B_TEMP(1,1,1) = BFX(BC%II-2,BC%JJ,BC%KK)
CALL GET_SCALAR_FACE_VALUE_NEW(U_TEMP,Z_TEMP,F_TEMP,B_TEMP,1,1,1,1,1,1,1,I_FLUX_LIMITER)
FX_H_S(BC%II-2,BC%JJ,BC%KK) = F_TEMP(1,1,1)
ENDIF
CASE( 2) OFF_WALL_SELECT_1
IF ((VV(BC%II,BC%JJ+1,BC%KK)>0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II,BC%JJ+1,BC%KK))%WALL_INDEX(+2)>0)) THEN
Z_TEMP(1,0:2,1) = (/RHO_H_S_P(BC%II,BC%JJ+1,BC%KK),RHO_H_S_P(BC%II,BC%JJ+1:BC%JJ+2,BC%KK)/)
U_TEMP(1,1,1) = VV(BC%II,BC%JJ+1,BC%KK)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,2,I_FLUX_LIMITER)
B_TEMP(1,1,1) = BFY(BC%II,BC%JJ+1,BC%KK)
CALL GET_SCALAR_FACE_VALUE_NEW(U_TEMP,Z_TEMP,F_TEMP,B_TEMP,1,1,1,1,1,1,2,I_FLUX_LIMITER)
FY_H_S(BC%II,BC%JJ+1,BC%KK) = F_TEMP(1,1,1)
ENDIF
CASE(-2) OFF_WALL_SELECT_1
IF ((VV(BC%II,BC%JJ-2,BC%KK)<0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II,BC%JJ-1,BC%KK))%WALL_INDEX(-2)>0)) THEN
Z_TEMP(1,1:3,1) = (/RHO_H_S_P(BC%II,BC%JJ-2:BC%JJ-1,BC%KK),RHO_H_S_P(BC%II,BC%JJ-1,BC%KK)/)
U_TEMP(1,1,1) = VV(BC%II,BC%JJ-2,BC%KK)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,2,I_FLUX_LIMITER)
B_TEMP(1,1,1) = BFY(BC%II,BC%JJ-2,BC%KK)
CALL GET_SCALAR_FACE_VALUE_NEW(U_TEMP,Z_TEMP,F_TEMP,B_TEMP,1,1,1,1,1,1,2,I_FLUX_LIMITER)
FY_H_S(BC%II,BC%JJ-2,BC%KK) = F_TEMP(1,1,1)
ENDIF
CASE( 3) OFF_WALL_SELECT_1
IF ((WW(BC%II,BC%JJ,BC%KK+1)>0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II,BC%JJ,BC%KK+1))%WALL_INDEX(+3)>0)) THEN
Z_TEMP(1,1,0:2) = (/RHO_H_S_P(BC%II,BC%JJ,BC%KK+1),RHO_H_S_P(BC%II,BC%JJ,BC%KK+1:BC%KK+2)/)
U_TEMP(1,1,1) = WW(BC%II,BC%JJ,BC%KK+1)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,3,I_FLUX_LIMITER)
B_TEMP(1,1,1) = BFZ(BC%II,BC%JJ,BC%KK+1)
CALL GET_SCALAR_FACE_VALUE_NEW(U_TEMP,Z_TEMP,F_TEMP,B_TEMP,1,1,1,1,1,1,3,I_FLUX_LIMITER)
FZ_H_S(BC%II,BC%JJ,BC%KK+1) = F_TEMP(1,1,1)
ENDIF
CASE(-3) OFF_WALL_SELECT_1
IF ((WW(BC%II,BC%JJ,BC%KK-2)<0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II,BC%JJ,BC%KK-1))%WALL_INDEX(-3)>0)) THEN
Z_TEMP(1,1,1:3) = (/RHO_H_S_P(BC%II,BC%JJ,BC%KK-2:BC%KK-1),RHO_H_S_P(BC%II,BC%JJ,BC%KK-1)/)
U_TEMP(1,1,1) = WW(BC%II,BC%JJ,BC%KK-2)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,3,I_FLUX_LIMITER)
B_TEMP(1,1,1) = BFZ(BC%II,BC%JJ,BC%KK-2)
CALL GET_SCALAR_FACE_VALUE_NEW(U_TEMP,Z_TEMP,F_TEMP,B_TEMP,1,1,1,1,1,1,3,I_FLUX_LIMITER)
FZ_H_S(BC%II,BC%JJ,BC%KK-2) = F_TEMP(1,1,1)
ENDIF
END SELECT OFF_WALL_SELECT_1
Expand Down Expand Up @@ -966,117 +976,7 @@ SUBROUTINE ENTHALPY_ADVECTION(U_DOT_DEL_RHO_H_S)
ENDDO
!$OMP END PARALLEL DO

END SUBROUTINE ENTHALPY_ADVECTION


SUBROUTINE SPECIES_ADVECTION_PART_1

USE MATH_FUNCTIONS, ONLY: GET_SCALAR_FACE_VALUE
USE PHYSICAL_FUNCTIONS, ONLY: GET_MOLECULAR_WEIGHT
REAL(EB), POINTER, DIMENSION(:,:,:) :: RHO_Z_P
REAL(EB), POINTER, DIMENSION(:,:,:,:) :: FX_ZZ,FY_ZZ,FZ_ZZ
REAL(EB), DIMENSION(0:3,0:3,0:3) :: U_TEMP,F_TEMP
REAL(EB), DIMENSION(-1:3,-1:3,-1:3) :: Z_TEMP
INTEGER :: I,J,K,IW,N
TYPE(WALL_TYPE), POINTER :: WC
TYPE(BOUNDARY_COORD_TYPE), POINTER :: BC
TYPE(BOUNDARY_PROP1_TYPE), POINTER :: B1

FX_ZZ=>SWORK1
FY_ZZ=>SWORK2
FZ_ZZ=>SWORK3
RHO_Z_P=>WORK_PAD

! Species face values

SPECIES_LOOP: DO N=1,N_TOTAL_SCALARS

!$OMP PARALLEL DO
DO K=-1,KBP1+1
DO J=-1,JBP1+1
DO I=-1,IBP1+1
RHO_Z_P(I,J,K) = RHOP(I,J,K)*ZZP(I,J,K,N)
ENDDO
ENDDO
ENDDO
!$OMP END PARALLEL DO

! Compute scalar face values

CALL GET_SCALAR_FACE_VALUE(UU,RHO_Z_P,FX_ZZ(:,:,:,N),0,IBAR,1,JBAR,1,KBAR,1,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_VALUE(VV,RHO_Z_P,FY_ZZ(:,:,:,N),1,IBAR,0,JBAR,1,KBAR,2,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_VALUE(WW,RHO_Z_P,FZ_ZZ(:,:,:,N),1,IBAR,1,JBAR,0,KBAR,3,I_FLUX_LIMITER)

!$OMP PARALLEL DO PRIVATE(IW,WC,BC,B1,U_TEMP,Z_TEMP,F_TEMP)
WALL_LOOP_2: DO IW=1,N_EXTERNAL_WALL_CELLS+N_INTERNAL_WALL_CELLS
WC=>WALL(IW)
IF (WC%BOUNDARY_TYPE==NULL_BOUNDARY) CYCLE WALL_LOOP_2
BC=>BOUNDARY_COORD(WC%BC_INDEX)
B1=>BOUNDARY_PROP1(WC%B1_INDEX)

! Overwrite first off-wall advective flux if flow is away from the wall and if the face is not also a wall cell

OFF_WALL_IF_2: IF (WC%BOUNDARY_TYPE/=INTERPOLATED_BOUNDARY .AND. WC%BOUNDARY_TYPE/=OPEN_BOUNDARY) THEN

OFF_WALL_SELECT_2: SELECT CASE(BC%IOR)
CASE( 1) OFF_WALL_SELECT_2
! ghost FX/UU(II+1)
! /// II /// II+1 | II+2 | ...
! ^ WALL_INDEX(II+1,+1)
IF ((UU(BC%II+1,BC%JJ,BC%KK)>0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II+1,BC%JJ,BC%KK))%WALL_INDEX(+1)>0)) THEN
Z_TEMP(0:2,1,1) = (/RHO_Z_P(BC%II+1,BC%JJ,BC%KK),RHO_Z_P(BC%II+1:BC%II+2,BC%JJ,BC%KK)/)
U_TEMP(1,1,1) = UU(BC%II+1,BC%JJ,BC%KK)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,1,I_FLUX_LIMITER)
FX_ZZ(BC%II+1,BC%JJ,BC%KK,N) = F_TEMP(1,1,1)
ENDIF
CASE(-1) OFF_WALL_SELECT_2
! FX/UU(II-2) ghost
! ... | II-2 | II-1 /// II ///
! ^ WALL_INDEX(II-1,-1)
IF ((UU(BC%II-2,BC%JJ,BC%KK)<0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II-1,BC%JJ,BC%KK))%WALL_INDEX(-1)>0)) THEN
Z_TEMP(1:3,1,1) = (/RHO_Z_P(BC%II-2:BC%II-1,BC%JJ,BC%KK),RHO_Z_P(BC%II-1,BC%JJ,BC%KK)/)
U_TEMP(1,1,1) = UU(BC%II-2,BC%JJ,BC%KK)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,1,I_FLUX_LIMITER)
FX_ZZ(BC%II-2,BC%JJ,BC%KK,N) = F_TEMP(1,1,1)
ENDIF
CASE( 2) OFF_WALL_SELECT_2
IF ((VV(BC%II,BC%JJ+1,BC%KK)>0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II,BC%JJ+1,BC%KK))%WALL_INDEX(+2)>0)) THEN
Z_TEMP(1,0:2,1) = (/RHO_Z_P(BC%II,BC%JJ+1,BC%KK),RHO_Z_P(BC%II,BC%JJ+1:BC%JJ+2,BC%KK)/)
U_TEMP(1,1,1) = VV(BC%II,BC%JJ+1,BC%KK)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,2,I_FLUX_LIMITER)
FY_ZZ(BC%II,BC%JJ+1,BC%KK,N) = F_TEMP(1,1,1)
ENDIF
CASE(-2) OFF_WALL_SELECT_2
IF ((VV(BC%II,BC%JJ-2,BC%KK)<0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II,BC%JJ-1,BC%KK))%WALL_INDEX(-2)>0)) THEN
Z_TEMP(1,1:3,1) = (/RHO_Z_P(BC%II,BC%JJ-2:BC%JJ-1,BC%KK),RHO_Z_P(BC%II,BC%JJ-1,BC%KK)/)
U_TEMP(1,1,1) = VV(BC%II,BC%JJ-2,BC%KK)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,2,I_FLUX_LIMITER)
FY_ZZ(BC%II,BC%JJ-2,BC%KK,N) = F_TEMP(1,1,1)
ENDIF
CASE( 3) OFF_WALL_SELECT_2
IF ((WW(BC%II,BC%JJ,BC%KK+1)>0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II,BC%JJ,BC%KK+1))%WALL_INDEX(+3)>0)) THEN
Z_TEMP(1,1,0:2) = (/RHO_Z_P(BC%II,BC%JJ,BC%KK+1),RHO_Z_P(BC%II,BC%JJ,BC%KK+1:BC%KK+2)/)
U_TEMP(1,1,1) = WW(BC%II,BC%JJ,BC%KK+1)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,3,I_FLUX_LIMITER)
FZ_ZZ(BC%II,BC%JJ,BC%KK+1,N) = F_TEMP(1,1,1)
ENDIF
CASE(-3) OFF_WALL_SELECT_2
IF ((WW(BC%II,BC%JJ,BC%KK-2)<0._EB) .AND. .NOT.(CELL(CELL_INDEX(BC%II,BC%JJ,BC%KK-1))%WALL_INDEX(-3)>0)) THEN
Z_TEMP(1,1,1:3) = (/RHO_Z_P(BC%II,BC%JJ,BC%KK-2:BC%KK-1),RHO_Z_P(BC%II,BC%JJ,BC%KK-1)/)
U_TEMP(1,1,1) = WW(BC%II,BC%JJ,BC%KK-2)
CALL GET_SCALAR_FACE_VALUE(U_TEMP,Z_TEMP,F_TEMP,1,1,1,1,1,1,3,I_FLUX_LIMITER)
FZ_ZZ(BC%II,BC%JJ,BC%KK-2,N) = F_TEMP(1,1,1)
ENDIF
END SELECT OFF_WALL_SELECT_2

ENDIF OFF_WALL_IF_2

ENDDO WALL_LOOP_2
!$OMP END PARALLEL DO

ENDDO SPECIES_LOOP

END SUBROUTINE SPECIES_ADVECTION_PART_1
END SUBROUTINE ENTHALPY_ADVECTION_NEW


SUBROUTINE SPECIES_ADVECTION_PART_1_NEW
Expand Down
29 changes: 17 additions & 12 deletions Source/init.f90
Original file line number Diff line number Diff line change
Expand Up @@ -595,18 +595,6 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM)

ALLOCATE(M%RSUM(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO) ; CALL ChkMemErr('INIT','RSUM',IZERO) ; M%RSUM=RSUM0

! Allocate scalar face values

ALLOCATE( M%FX(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FX',IZERO) ; M%FX=0._EB
ALLOCATE( M%FY(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FY',IZERO) ; M%FY=0._EB
ALLOCATE( M%FZ(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FZ',IZERO) ; M%FZ=0._EB

IF (FLUX_LIMITER_SINGLE_COEF) THEN
ALLOCATE( M%BFX(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO); CALL ChkMemErr('INIT','BFX',IZERO) ; M%BFX=0._EB
ALLOCATE( M%BFY(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO); CALL ChkMemErr('INIT','BFY',IZERO) ; M%BFY=0._EB
ALLOCATE( M%BFZ(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO); CALL ChkMemErr('INIT','BFZ',IZERO) ; M%BFZ=0._EB
ENDIF

! Allocate storage for scalar total fluxes

IF (STORE_SPECIES_FLUX) THEN
Expand Down Expand Up @@ -831,6 +819,23 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM)
M%ZZS(:,:,:,N) = INITIAL_UNMIXED_FRACTION
ENDDO

! Allocate scalar face values

ALLOCATE( M%FX(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FX',IZERO) ; M%FX=0._EB
ALLOCATE( M%FY(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FY',IZERO) ; M%FY=0._EB
ALLOCATE( M%FZ(0:IBP1,0:JBP1,0:KBP1,1:N_TOTAL_SCALARS),STAT=IZERO); CALL ChkMemErr('INIT','FZ',IZERO) ; M%FZ=0._EB
DO N=1,N_TRACKED_SPECIES
M%FX(:,:,:,N)=RHOA*SPECIES_MIXTURE(N)%ZZ0
M%FY(:,:,:,N)=RHOA*SPECIES_MIXTURE(N)%ZZ0
M%FZ(:,:,:,N)=RHOA*SPECIES_MIXTURE(N)%ZZ0
ENDDO

! Allocate flux limiter coefficients

ALLOCATE( M%BFX(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO); CALL ChkMemErr('INIT','BFX',IZERO) ; M%BFX=0._EB
ALLOCATE( M%BFY(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO); CALL ChkMemErr('INIT','BFY',IZERO) ; M%BFY=0._EB
ALLOCATE( M%BFZ(0:IBP1,0:JBP1,0:KBP1),STAT=IZERO); CALL ChkMemErr('INIT','BFZ',IZERO) ; M%BFZ=0._EB

! Allocate and Initialize Mesh-Dependent Radiation Arrays

M%QR = 0._EB
Expand Down
12 changes: 2 additions & 10 deletions Source/main.f90
Original file line number Diff line number Diff line change
Expand Up @@ -626,11 +626,7 @@ PROGRAM FDS
COMPUTE_FINITE_DIFFERENCES_1: DO NM=LOWER_MESH_INDEX,UPPER_MESH_INDEX
CALL INSERT_ALL_PARTICLES(T,NM)
IF (.NOT.SOLID_PHASE_ONLY .AND. .NOT.FREEZE_VELOCITY) CALL COMPUTE_VISCOSITY(NM,APPLY_TO_ESTIMATED_VARIABLES=.FALSE.)
IF (FLUX_LIMITER_SINGLE_COEF) THEN
CALL MASS_FINITE_DIFFERENCES_NEW(NM)
ELSE
CALL MASS_FINITE_DIFFERENCES(NM)
ENDIF
CALL MASS_FINITE_DIFFERENCES_NEW(NM)
ENDDO COMPUTE_FINITE_DIFFERENCES_1

! Estimate quantities at next time step, and decrease/increase time step if necessary based on CFL condition
Expand Down Expand Up @@ -811,11 +807,7 @@ PROGRAM FDS

COMPUTE_FINITE_DIFFERENCES_2: DO NM=LOWER_MESH_INDEX,UPPER_MESH_INDEX
IF (.NOT.SOLID_PHASE_ONLY .AND. .NOT.FREEZE_VELOCITY) CALL COMPUTE_VISCOSITY(NM,APPLY_TO_ESTIMATED_VARIABLES=.TRUE.)
IF (FLUX_LIMITER_SINGLE_COEF) THEN
CALL MASS_FINITE_DIFFERENCES_NEW(NM)
ELSE
CALL MASS_FINITE_DIFFERENCES(NM)
ENDIF
CALL MASS_FINITE_DIFFERENCES_NEW(NM)
CALL DENSITY(T,DT,NM)
IF (LEVEL_SET_MODE>0) CALL LEVEL_SET_FIRESPREAD(T,DT,NM)
ENDDO COMPUTE_FINITE_DIFFERENCES_2
Expand Down
Loading
Loading