Skip to content

Commit

Permalink
FDS Source: Revert change to divg.f90 made in commit cccc20d
Browse files Browse the repository at this point in the history
  • Loading branch information
mcgratta committed Jun 16, 2022
1 parent 4bc9f00 commit 5984999
Show file tree
Hide file tree
Showing 2 changed files with 76 additions and 56 deletions.
128 changes: 74 additions & 54 deletions Source/divg.f90
Original file line number Diff line number Diff line change
Expand Up @@ -821,7 +821,7 @@ SUBROUTINE DIVERGENCE_PART_1(T,DT,NM)
SUBROUTINE ENTHALPY_ADVECTION

REAL(EB), POINTER, DIMENSION(:,:,:) :: FX_H_S,FY_H_S,FZ_H_S
REAL(EB) :: UN_P,TMP_F_GAS,DU_P,DU_M,DV_P,DV_M,DW_P,DW_M
REAL(EB) :: UN,UN_P,TMP_F_GAS,DU_P,DU_M,DV_P,DV_M,DW_P,DW_M,DU
INTEGER :: IC

RHO_H_S_P=>WORK1 ; RHO_H_S_P = 0._EB
Expand Down Expand Up @@ -854,11 +854,8 @@ SUBROUTINE ENTHALPY_ADVECTION
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)

!$OMP PARALLEL PRIVATE(ZZ_GET)

ALLOCATE(ZZ_GET(1:N_TRACKED_SPECIES))

!$OMP DO PRIVATE(IW,WC,BC,ONE_D,II,JJ,KK,IIG,JJG,KKG,IOR,UN_P,TMP_F_GAS,H_S,Z_TEMP,U_TEMP,F_TEMP,IC)
WALL_LOOP: DO IW=1,N_EXTERNAL_WALL_CELLS+N_INTERNAL_WALL_CELLS

WC=>WALL(IW)
Expand Down Expand Up @@ -948,50 +945,63 @@ SUBROUTINE ENTHALPY_ADVECTION

ENDIF

IC = CELL_INDEX(II,JJ,KK)
IF ((SOLID(IC) .AND. .NOT.EXTERIOR(IC)) .OR. EXTERIOR(IC)) THEN
SELECT CASE(IOR)
CASE( 1); FX_H_S(II,JJ,KK) = ONE_D%RHO_F*H_S
CASE(-1); FX_H_S(II-1,JJ,KK) = ONE_D%RHO_F*H_S
CASE( 2); FY_H_S(II,JJ,KK) = ONE_D%RHO_F*H_S
CASE(-2); FY_H_S(II,JJ-1,KK) = ONE_D%RHO_F*H_S
CASE( 3); FZ_H_S(II,JJ,KK) = ONE_D%RHO_F*H_S
CASE(-3); FZ_H_S(II,JJ,KK-1) = ONE_D%RHO_F*H_S
END SELECT
ENDIF
SELECT CASE(WC%BOUNDARY_TYPE)
CASE DEFAULT
IOR_SELECT: SELECT CASE(IOR)
CASE( 1); UN = UU(II,JJ,KK)
CASE(-1); UN = UU(II-1,JJ,KK)
CASE( 2); UN = VV(II,JJ,KK)
CASE(-2); UN = VV(II,JJ-1,KK)
CASE( 3); UN = WW(II,JJ,KK)
CASE(-3); UN = WW(II,JJ,KK-1)
END SELECT IOR_SELECT
CASE(SOLID_BOUNDARY)
IF (PREDICTOR) UN = -SIGN(1._EB,REAL(IOR,EB))*ONE_D%U_NORMAL_S
IF (CORRECTOR) UN = -SIGN(1._EB,REAL(IOR,EB))*ONE_D%U_NORMAL
CASE(INTERPOLATED_BOUNDARY)
UN = UVW_SAVE(IW)
END SELECT

DU = (ONE_D%RHO_F*H_S - RHO_H_S_P(IIG,JJG,KKG))*UN
U_DOT_DEL_RHO_H_S(IIG,JJG,KKG) = U_DOT_DEL_RHO_H_S(IIG,JJG,KKG) - SIGN(1._EB,REAL(IOR,EB))*DU*ONE_D%RDN

ENDDO WALL_LOOP
!$OMP END DO

DEALLOCATE(ZZ_GET)

! FDS Tech Guide (B.12-B.14)

!$OMP DO PRIVATE(DU_P,DU_M,DV_P,DV_M,DW_P,DW_M)
!$OMP PARALLEL DO PRIVATE(IC,DU_P,DU_M,DV_P,DV_M,DW_P,DW_M)
DO K=1,KBAR
DO J=1,JBAR
DO I=1,IBAR
DU_P = (FX_H_S(I,J,K) - RHO_H_S_P(I,J,K))*UU(I,J,K)
DU_M = (FX_H_S(I-1,J,K) - RHO_H_S_P(I,J,K))*UU(I-1,J,K)
DV_P = (FY_H_S(I,J,K) - RHO_H_S_P(I,J,K))*VV(I,J,K)
DV_M = (FY_H_S(I,J-1,K) - RHO_H_S_P(I,J,K))*VV(I,J-1,K)
DW_P = (FZ_H_S(I,J,K) - RHO_H_S_P(I,J,K))*WW(I,J,K)
DW_M = (FZ_H_S(I,J,K-1) - RHO_H_S_P(I,J,K))*WW(I,J,K-1)
U_DOT_DEL_RHO_H_S(I,J,K) = (DU_P-DU_M)*RDX(I) + (DV_P-DV_M)*RDY(J) + (DW_P-DW_M)*RDZ(K)
IC = CELL_INDEX(I,J,K)
IF (SOLID(IC)) CYCLE
DU_P = 0._EB
DU_M = 0._EB
DV_P = 0._EB
DV_M = 0._EB
DW_P = 0._EB
DW_M = 0._EB
IF (WALL_INDEX(IC, 1)==0) DU_P = (FX_H_S(I,J,K) - RHO_H_S_P(I,J,K))*UU(I,J,K)
IF (WALL_INDEX(IC,-1)==0) DU_M = (FX_H_S(I-1,J,K) - RHO_H_S_P(I,J,K))*UU(I-1,J,K)
IF (WALL_INDEX(IC, 2)==0) DV_P = (FY_H_S(I,J,K) - RHO_H_S_P(I,J,K))*VV(I,J,K)
IF (WALL_INDEX(IC,-2)==0) DV_M = (FY_H_S(I,J-1,K) - RHO_H_S_P(I,J,K))*VV(I,J-1,K)
IF (WALL_INDEX(IC, 3)==0) DW_P = (FZ_H_S(I,J,K) - RHO_H_S_P(I,J,K))*WW(I,J,K)
IF (WALL_INDEX(IC,-3)==0) DW_M = (FZ_H_S(I,J,K-1) - RHO_H_S_P(I,J,K))*WW(I,J,K-1)
U_DOT_DEL_RHO_H_S(I,J,K) = U_DOT_DEL_RHO_H_S(I,J,K) + (DU_P-DU_M)*RDX(I) + (DV_P-DV_M)*RDY(J) + (DW_P-DW_M)*RDZ(K)
ENDDO
ENDDO
ENDDO
!$OMP END DO

!$OMP END PARALLEL
!$OMP END PARALLEL DO

END SUBROUTINE ENTHALPY_ADVECTION


SUBROUTINE SPECIES_ADVECTION

REAL(EB), POINTER, DIMENSION(:,:,:) :: FX_ZZ,FY_ZZ,FZ_ZZ
REAL(EB) :: DU_P,DU_M,DV_P,DV_M,DW_P,DW_M
REAL(EB) :: UN,DU_P,DU_M,DV_P,DV_M,DW_P,DW_M,DU
INTEGER :: IC

FX_ZZ=>WORK2 ; FX_ZZ = 0._EB
Expand All @@ -1016,9 +1026,6 @@ SUBROUTINE SPECIES_ADVECTION
CALL GET_SCALAR_FACE_VALUE(VV,RHO_Z_P,FY_ZZ,1,IBAR,1,JBM1,1,KBAR,2,I_FLUX_LIMITER)
CALL GET_SCALAR_FACE_VALUE(WW,RHO_Z_P,FZ_ZZ,1,IBAR,1,JBAR,1,KBM1,3,I_FLUX_LIMITER)

!$OMP PARALLEL

!$OMP DO SCHEDULE(DYNAMIC) PRIVATE(IW,WC,ONE_D,BC,II,JJ,KK,IIG,JJG,KKG,IOR,U_TEMP,Z_TEMP,F_TEMP,IC)
WALL_LOOP: DO IW=1,N_EXTERNAL_WALL_CELLS+N_INTERNAL_WALL_CELLS

WC=>WALL(IW)
Expand Down Expand Up @@ -1092,38 +1099,51 @@ SUBROUTINE SPECIES_ADVECTION

! Correct flux terms at the boundary

IC = CELL_INDEX(II,JJ,KK)
IF ((SOLID(IC) .AND. .NOT.EXTERIOR(IC)) .OR. EXTERIOR(IC)) THEN
SELECT CASE(IOR)
CASE( 1); FX_ZZ(II,JJ,KK) = ONE_D%RHO_F*ONE_D%ZZ_F(N)
CASE(-1); FX_ZZ(II-1,JJ,KK) = ONE_D%RHO_F*ONE_D%ZZ_F(N)
CASE( 2); FY_ZZ(II,JJ,KK) = ONE_D%RHO_F*ONE_D%ZZ_F(N)
CASE(-2); FY_ZZ(II,JJ-1,KK) = ONE_D%RHO_F*ONE_D%ZZ_F(N)
CASE( 3); FZ_ZZ(II,JJ,KK) = ONE_D%RHO_F*ONE_D%ZZ_F(N)
CASE(-3); FZ_ZZ(II,JJ,KK-1) = ONE_D%RHO_F*ONE_D%ZZ_F(N)
END SELECT
ENDIF
SELECT CASE(WC%BOUNDARY_TYPE)
CASE DEFAULT
IOR_SELECT: SELECT CASE(IOR)
CASE( 1); UN = UU(II,JJ,KK)
CASE(-1); UN = UU(II-1,JJ,KK)
CASE( 2); UN = VV(II,JJ,KK)
CASE(-2); UN = VV(II,JJ-1,KK)
CASE( 3); UN = WW(II,JJ,KK)
CASE(-3); UN = WW(II,JJ,KK-1)
END SELECT IOR_SELECT
CASE(SOLID_BOUNDARY)
IF (PREDICTOR) UN = -SIGN(1._EB,REAL(IOR,EB))*ONE_D%U_NORMAL_S
IF (CORRECTOR) UN = -SIGN(1._EB,REAL(IOR,EB))*ONE_D%U_NORMAL
CASE(INTERPOLATED_BOUNDARY)
UN = UVW_SAVE(IW)
END SELECT

DU = (ONE_D%RHO_F*ONE_D%ZZ_F(N) - RHO_Z_P(IIG,JJG,KKG))*UN
U_DOT_DEL_RHO_Z(IIG,JJG,KKG) = U_DOT_DEL_RHO_Z(IIG,JJG,KKG) - SIGN(1._EB,REAL(IOR,EB))*DU*ONE_D%RDN

ENDDO WALL_LOOP
!$OMP END DO

!$OMP DO PRIVATE(DU_P,DU_M,DV_P,DV_M,DW_P,DW_M) SCHEDULE(STATIC)
!$OMP PARALLEL DO PRIVATE(IC,DU_P,DU_M,DV_P,DV_M,DW_P,DW_M)
DO K=1,KBAR
DO J=1,JBAR
DO I=1,IBAR
DU_P = (FX_ZZ(I,J,K) - RHO_Z_P(I,J,K))*UU(I,J,K)
DU_M = (FX_ZZ(I-1,J,K) - RHO_Z_P(I,J,K))*UU(I-1,J,K)
DV_P = (FY_ZZ(I,J,K) - RHO_Z_P(I,J,K))*VV(I,J,K)
DV_M = (FY_ZZ(I,J-1,K) - RHO_Z_P(I,J,K))*VV(I,J-1,K)
DW_P = (FZ_ZZ(I,J,K) - RHO_Z_P(I,J,K))*WW(I,J,K)
DW_M = (FZ_ZZ(I,J,K-1) - RHO_Z_P(I,J,K))*WW(I,J,K-1)
U_DOT_DEL_RHO_Z(I,J,K) = (DU_P-DU_M)*RDX(I) + (DV_P-DV_M)*RDY(J) + (DW_P-DW_M)*RDZ(K)
IC = CELL_INDEX(I,J,K)
IF (SOLID(IC)) CYCLE
DU_P = 0._EB
DU_M = 0._EB
DV_P = 0._EB
DV_M = 0._EB
DW_P = 0._EB
DW_M = 0._EB
IF (WALL_INDEX(IC, 1)==0) DU_P = (FX_ZZ(I,J,K) - RHO_Z_P(I,J,K))*UU(I,J,K)
IF (WALL_INDEX(IC,-1)==0) DU_M = (FX_ZZ(I-1,J,K) - RHO_Z_P(I,J,K))*UU(I-1,J,K)
IF (WALL_INDEX(IC, 2)==0) DV_P = (FY_ZZ(I,J,K) - RHO_Z_P(I,J,K))*VV(I,J,K)
IF (WALL_INDEX(IC,-2)==0) DV_M = (FY_ZZ(I,J-1,K) - RHO_Z_P(I,J,K))*VV(I,J-1,K)
IF (WALL_INDEX(IC,+3)==0) DW_P = (FZ_ZZ(I,J,K) - RHO_Z_P(I,J,K))*WW(I,J,K)
IF (WALL_INDEX(IC,-3)==0) DW_M = (FZ_ZZ(I,J,K-1) - RHO_Z_P(I,J,K))*WW(I,J,K-1)
U_DOT_DEL_RHO_Z(I,J,K) = U_DOT_DEL_RHO_Z(I,J,K) + (DU_P-DU_M)*RDX(I) + (DV_P-DV_M)*RDY(J) + (DW_P-DW_M)*RDZ(K)
ENDDO
ENDDO
ENDDO
!$OMP END DO

!$OMP END PARALLEL
!$OMP END PARALLEL DO

END SUBROUTINE SPECIES_ADVECTION

Expand Down
4 changes: 2 additions & 2 deletions Verification/WUI/char_oxidation_1.fds
Original file line number Diff line number Diff line change
Expand Up @@ -65,13 +65,13 @@

&INIT ID = 'foliage part',
PART_ID = 'foliage part'
XB = -0.3, 0.1,-0.2, 0.2,0.4,0.6
XB = -0.3, 0.1,-0.2, 0.2,0.2,0.4
N_PARTICLES_PER_CELL = 1
CELL_CENTERED = .TRUE.
MASS_PER_VOLUME = 1.66
DRY = T /

&SURF ID='radiant panel', TMP_FRONT=1800, COLOR='RED' /
&SURF ID='radiant panel', TMP_FRONT=1700, COLOR='RED' /

&VENT MB='XMIN',SURF_ID='radiant panel' /
&VENT MB='XMAX',SURF_ID='OPEN' /
Expand Down

0 comments on commit 5984999

Please sign in to comment.