diff --git a/Source/geom.f90 b/Source/geom.f90 index 04187148d15..919f11e4e44 100644 --- a/Source/geom.f90 +++ b/Source/geom.f90 @@ -293,6 +293,15 @@ MODULE COMPLEX_GEOMETRY ! Areas per SURF and GEOM: REAL(EB), ALLOCATABLE, DIMENSION(:,:) :: FDS_AREA_GEOM +! Cut-cell state and CFACE origin parameter: +INTEGER, PARAMETER :: NOT_BLOCKED = 0 +INTEGER, PARAMETER :: BLOCKED_SMALL_CELL = 1 +INTEGER, PARAMETER :: BLOCKED_SPLIT_CELL = 2 +INTEGER, PARAMETER :: BLOCKED_REFI_INTER = 3 +INTEGER, PARAMETER :: BLOCKED_CAVITY_CELL= 4 +INTEGER, PARAMETER :: BLOCKED_UNLINK_CELL= 5 + + ! End Variable declaration for CC_IBM. !! --------------------------------------------------------------------------------- @@ -747,7 +756,7 @@ SUBROUTINE SET_CUTCELLS_3D TYPE(WALL_TYPE), POINTER :: WC TYPE(BOUNDARY_COORD_TYPE), POINTER :: BC INTEGER, ALLOCATABLE, DIMENSION(:,:) :: CCELEM,FACE_LIST -LOGICAL, ALLOCATABLE, DIMENSION(:) :: NOADVANCE +INTEGER, ALLOCATABLE, DIMENSION(:) :: NOADVANCE REAL(EB),ALLOCATABLE, DIMENSION(:,:) :: XYZCEN REAL(EB),ALLOCATABLE, DIMENSION(:) :: VOLUME INTEGER, PARAMETER :: ADDI(LOW_IND:HIGH_IND,IAXIS:KAXIS) = RESHAPE((/-1,0, 0,0, 0,0/),(/2,3/)) @@ -1394,6 +1403,10 @@ SUBROUTINE SET_CUTCELLS_3D CASE(-KAXIS,KAXIS); VOLUME(1:CC%NCELL) = CC%XYZCEN(KAXIS,1:CC%NCELL) CASE DEFAULT; VOLUME(1:CC%NCELL) = CC%VOLUME(1:CC%NCELL) END SELECT + ! IF(NORM2(CELL_BLOCK_NVEC)>TWO_EPSILON_EB) THEN ! Block in the NVEC direction: + ! DO J=1,CC%NCELL; VOLUME(J) = DOT_PRODUCT(CELL_BLOCK_NVEC,CC%XYZCEN(IAXIS:KAXIS,J)); ENDDO + ! ENDIF + DO J=1,CC%NCELL; VOLUME(J) = VOLUME(J)*(1._EB+REAL(J-1,EB)*GEOMEPS); ENDDO SELECT CASE(CELL_BLOCK_IOR) CASE(-IAXIS); I=MAXLOC(VOLUME(1:CC%NCELL),DIM=1) @@ -1407,7 +1420,7 @@ SUBROUTINE SET_CUTCELLS_3D DEALLOCATE(VOLUME) NCELL_LOOP_2 : DO J=1,CC%NCELL IF(J==I) CYCLE NCELL_LOOP_2 - CC%NOADVANCE(J)=.TRUE. + CC%NOADVANCE(J) = BLOCKED_SPLIT_CELL ENDDO NCELL_LOOP_2 ENDDO ENDIF @@ -1425,7 +1438,7 @@ SUBROUTINE SET_CUTCELLS_3D ICC = MESHES(NM)%CCVAR(I,J,K,CC_IDCC); IF(ICC<1) CYCLE CC =>MESHES(NM)%CUT_CELL(ICC) DO JCC=1,CC%NCELL - IF(.NOT.CC%NOADVANCE(JCC)) CYCLE + IF(CC%NOADVANCE(JCC)<1) CYCLE DO IFC=2,CC%CCELEM(1,JCC)+1 IFACE=CC%CCELEM(IFC,JCC) IF(CC%FACE_LIST(1,IFACE)/=CC_FTYPE_CFINB) CYCLE @@ -1460,7 +1473,7 @@ SUBROUTINE SET_CUTCELLS_3D END SELECT ENDDO IF(SUM_FACE>1 .OR. SUM_CCELL>0) CYCLE - CC%NOADVANCE(J)=.TRUE. + CC%NOADVANCE(J)=BLOCKED_CAVITY_CELL K=K+1 ENDDO ENDDO @@ -2236,7 +2249,7 @@ SUBROUTINE TAG_CC_BLOCKING_REFINEMENT M => MESHES(NM) ! Set all fine side cut-cells in cells next to external boundaries which have SOLID coarse mesh faces - ! to CC%NOADVANCE(J)=.TRUE. and block them. + ! to CC%NOADVANCE(J)=BLOCKED_REFI_INTER and block them. EXT_WALL_LOOP_1 : DO IW=1,M%N_EXTERNAL_WALL_CELLS WC=>WALL(IW) EWC=>EXTERNAL_WALL(IW) @@ -2372,7 +2385,7 @@ SUBROUTINE TAG_CC_BLOCKING_REFINEMENT M => MESHES(NM) ! Set all fine side cut-cells in cells next to external boundaries which have SOLID coarse mesh faces - ! to CC%NOADVANCE(J)=.TRUE. and block them. + ! to CC%NOADVANCE(J)=BLOCKED_REFI_INTER and block them. EXT_WALL_LOOP : DO IW=1,M%N_EXTERNAL_WALL_CELLS WC=>WALL(IW); IF (WC%BOUNDARY_TYPE/=INTERPOLATED_BOUNDARY) CYCLE EXT_WALL_LOOP EWC=>EXTERNAL_WALL(IW) @@ -2408,7 +2421,7 @@ SUBROUTINE TAG_CC_BLOCKING_REFINEMENT ALLOCATE(FACE_LIST(1:CC_NPARAM_CCFACE,1:NFACE_CELL)); FACE_LIST = CC_UNDEFINED ALLOCATE(VOLUME(1:NCELL)); VOLUME(1)=M%DX(II)*M%DY(JJ)*M%DZ(KK) ALLOCATE(XYZCEN(IAXIS:KAXIS,1:NCELL)); XYZCEN(IAXIS:KAXIS,1) = (/ M%XC(II),M%YC(JJ),M%ZC(KK) /) - ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = .TRUE. + ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = BLOCKED_REFI_INTER ! Add one by one regular and gas cut faces: CT = 1; CCELEM(1,1) = 0 DO AX=IAXIS,KAXIS @@ -2442,7 +2455,7 @@ SUBROUTINE TAG_CC_BLOCKING_REFINEMENT IFACE = CC%CCELEM(IFC,JCC) IF( (CC%FACE_LIST(1,IFACE)==CC_FTYPE_RCGAS .OR. CC%FACE_LIST(1,IFACE)==CC_FTYPE_CFGAS) .AND. & CC%FACE_LIST(2,IFACE)==LOHIF .AND. CC%FACE_LIST(3,IFACE)==X1AXIS ) THEN - CC%NOADVANCE(JCC)=.TRUE. + CC%NOADVANCE(JCC)=BLOCKED_REFI_INTER CYCLE JCC_LOOP_1 ENDIF ENDDO @@ -2480,7 +2493,7 @@ SUBROUTINE TAG_CC_BLOCKING_REFINEMENT ALLOCATE(FACE_LIST(1:CC_NPARAM_CCFACE,1:NFACE_CELL)); FACE_LIST = CC_UNDEFINED ALLOCATE(VOLUME(1:NCELL)); VOLUME(1)=M2%DX(IOGC)*M2%DY(JOGC)*M2%DZ(KOGC) ALLOCATE(XYZCEN(IAXIS:KAXIS,1:NCELL)); XYZCEN(IAXIS:KAXIS,1) = (/ M2%XC(IOGC),M2%YC(JOGC),M2%ZC(KOGC) /) - ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = .TRUE. + ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = BLOCKED_REFI_INTER ! Add one by one regular and gas cut faces: CT = 1; CCELEM(1,1) = 0 DO AX=IAXIS,KAXIS @@ -2515,7 +2528,7 @@ SUBROUTINE TAG_CC_BLOCKING_REFINEMENT IFACE = CC%CCELEM(IFC,JCC) IF( (CC%FACE_LIST(1,IFACE)==CC_FTYPE_RCGAS .OR. CC%FACE_LIST(1,IFACE)==CC_FTYPE_CFGAS) .AND. & CC%FACE_LIST(2,IFACE)==LOHIF .AND. CC%FACE_LIST(3,IFACE)==X1AXIS ) THEN - CC%NOADVANCE(JCC)=.TRUE. + CC%NOADVANCE(JCC)=BLOCKED_REFI_INTER CYCLE JCC_LOOP_3 ENDIF ENDDO @@ -2557,7 +2570,7 @@ SUBROUTINE TAG_BLOCK_CELL(NM,II1,JJ1,KK1,NOM,IIO1,JJO1,KKO1,FINE_CELL) ALLOCATE(FACE_LIST(1:CC_NPARAM_CCFACE,1:NFACE_CELL)); FACE_LIST = CC_UNDEFINED ALLOCATE(VOLUME(1:NCELL)); VOLUME(1)=M%DX(II1)*M%DY(JJ1)*M%DZ(KK1) ALLOCATE(XYZCEN(IAXIS:KAXIS,1:NCELL)); XYZCEN(IAXIS:KAXIS,1) = (/ M%XC(II1),M%YC(JJ1),M%ZC(KK1) /) - ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = .FALSE. + ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = BLOCKED_REFI_INTER ! Add one by one regular and gas cut faces: CT = 1; CCELEM(1,1) = 0 DO AX=IAXIS,KAXIS @@ -2588,7 +2601,8 @@ SUBROUTINE TAG_BLOCK_CELL(NM,II1,JJ1,KK1,NOM,IIO1,JJO1,KKO1,FINE_CELL) ENDIF ! Here Test if cut-cells in II,KK,KK are blocked or not in IIO,JJO,KKO: IF(ICC>0) THEN - IF(M2%CCVAR(IIO1,JJO1,KKO1,CC_CGSC)==CC_SOLID) THEN; M%CUT_CELL(ICC)%NOADVANCE(1:M%CUT_CELL(ICC)%NCELL) = .TRUE. + IF(M2%CCVAR(IIO1,JJO1,KKO1,CC_CGSC)==CC_SOLID) THEN + M%CUT_CELL(ICC)%NOADVANCE(1:M%CUT_CELL(ICC)%NCELL) = BLOCKED_REFI_INTER ELSE; CALL TEST_CC_FOR_BLOCKING(NM,ICC,NOM,IIO1,JJO1,KKO1,ICC2) ENDIF ENDIF @@ -2609,7 +2623,7 @@ SUBROUTINE TAG_BLOCK_CELL(NM,II1,JJ1,KK1,NOM,IIO1,JJO1,KKO1,FINE_CELL) ALLOCATE(FACE_LIST(1:CC_NPARAM_CCFACE,1:NFACE_CELL)); FACE_LIST = CC_UNDEFINED ALLOCATE(VOLUME(1:NCELL)); VOLUME(1)=M2%DX(IIO1)*M2%DY(JJO1)*M2%DZ(KKO1) ALLOCATE(XYZCEN(IAXIS:KAXIS,1:NCELL)); XYZCEN(IAXIS:KAXIS,1) = (/ M2%XC(IIO1),M2%YC(JJO1),M2%ZC(KKO1) /) - ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = .FALSE. + ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = NOT_BLOCKED ! Add one by one regular and gas cut faces: CT = 1; CCELEM(1,1) = 0 DO AX=IAXIS,KAXIS @@ -2907,7 +2921,7 @@ SUBROUTINE TEST_CC_FOR_BLOCKING(NM,ICC,NOM,IIO1,JJO1,KKO1,ICC2) ! ENDIF ENDDO INBFC_LOC_LOOP ! Here set no ADVANCE if BLOCK_CELL=T: - IF(BLOCK_CELL) CC%NOADVANCE(JCC)=.TRUE. + IF(BLOCK_CELL) CC%NOADVANCE(JCC) = BLOCKED_REFI_INTER ENDDO JCC_LOOP ! IF(NM==1 .AND. ICC<30) CLOSE(LU_CCELL) @@ -4145,6 +4159,7 @@ SUBROUTINE GET_EXT_INB_CUTFACES_TO_CFACE SURF_INDEX = CF%SURF_INDEX(IFACE) CALL INIT_CFACE_CELL(NM,ICF,IFACE,CFACE_INDEX_LOCAL,SURF_INDEX,INTEGER_ONE,IS_INB=.TRUE.) ENDDO + IF(ALLOCATED(CF%CFACE_ORIGIN)) DEALLOCATE(CF%CFACE_ORIGIN) ENDDO MESHES(NM)%N_INTERNAL_CFACE_CELLS = CFACE_INDEX_LOCAL - MESHES(NM)%INTERNAL_CFACE_CELLS_LB ENDDO MESH_LOOP_1 @@ -4863,7 +4878,7 @@ SUBROUTINE BLOCK_SMALL_UNLINKED_CUTCELLS(NM,NBLKCELLS) ICC = M%CCVAR(I,J,K,CC_IDCC) IF(ICC<1) CYCLE ! No Cut-cell. JCC_LOOP : DO JCC=1,M%CUT_CELL(ICC)%NCELL - IF(.NOT.M%CUT_CELL(ICC)%NOADVANCE(JCC) .OR. M%CUT_CELL(ICC)%IJK_LINK(1,JCC)/=CC_SOLID) CYCLE JCC_LOOP + IF(M%CUT_CELL(ICC)%NOADVANCE(JCC)<1 .OR. M%CUT_CELL(ICC)%IJK_LINK(1,JCC)/=CC_SOLID) CYCLE JCC_LOOP NBLKCELLS = NBLKCELLS + 1 CALL BLOCK_CUT_CELL(NM,ICC,JCC,1) ENDDO JCC_LOOP @@ -4880,7 +4895,7 @@ SUBROUTINE BLOCK_SMALL_UNLINKED_CUTCELLS(NM,NBLKCELLS) IF(ICC<1) CYCLE ! No Cut-cell. NCELL = M%CUT_CELL(ICC)%NCELL JCC_LOOP_2 : DO JCC=1,NCELL - IF(.NOT.M%CUT_CELL(ICC)%NOADVANCE(JCC) .OR. M%CUT_CELL(ICC)%IJK_LINK(1,JCC)/=CC_SOLID) CYCLE JCC_LOOP_2 + IF(M%CUT_CELL(ICC)%NOADVANCE(JCC)<1 .OR. M%CUT_CELL(ICC)%IJK_LINK(1,JCC)/=CC_SOLID) CYCLE JCC_LOOP_2 CALL BLOCK_CUT_CELL(NM,ICC,JCC,2) ENDDO JCC_LOOP_2 ENDDO @@ -4895,7 +4910,7 @@ SUBROUTINE BLOCK_SMALL_UNLINKED_CUTCELLS(NM,NBLKCELLS) IF(ICC<1) CYCLE ! No Cut-cell. NCELL = M%CUT_CELL(ICC)%NCELL JCC_LOOP_3 : DO JCC=NCELL,1,-1 - IF(.NOT.M%CUT_CELL(ICC)%NOADVANCE(JCC) .OR. M%CUT_CELL(ICC)%IJK_LINK(1,JCC)/=CC_SOLID) CYCLE JCC_LOOP_3 + IF(M%CUT_CELL(ICC)%NOADVANCE(JCC)<1 .OR. M%CUT_CELL(ICC)%IJK_LINK(1,JCC)/=CC_SOLID) CYCLE JCC_LOOP_3 CALL BLOCK_CUT_CELL(NM,ICC,JCC,3) ENDDO JCC_LOOP_3 ENDDO @@ -5080,7 +5095,7 @@ SUBROUTINE GET_REMAINING_CUTCELLS(NM) INTEGER :: I,J,K,CT,X1AXIS,SIDE,ICC,JCC,IFACE,ICF,JCF,ICFC,ICFINB,NCFACE_CUTCELL,NCELL,NFACE_CELL INTEGER :: NCC_MESH,NGC_MESH,NCELL_IN,NCELL_GC,COUNT_CC,COUNT_GC INTEGER, ALLOCATABLE, DIMENSION(:,:) :: CCELEM,FACE_LIST -LOGICAL, ALLOCATABLE, DIMENSION(:) :: NOADVANCE +INTEGER, ALLOCATABLE, DIMENSION(:) :: NOADVANCE REAL(EB),ALLOCATABLE, DIMENSION(:,:) :: XYZCEN REAL(EB),ALLOCATABLE, DIMENSION(:) :: VOLUME INTEGER, PARAMETER :: ADDI(LOW_IND:HIGH_IND,IAXIS:KAXIS) = RESHAPE((/-1,0, 0,0, 0,0/),(/2,3/)) @@ -5190,7 +5205,7 @@ SUBROUTINE GET_REMAINING_CUTCELLS(NM) ALLOCATE(FACE_LIST(1:CC_NPARAM_CCFACE,1:NFACE_CELL)); FACE_LIST = CC_UNDEFINED ALLOCATE(VOLUME(1:NCELL)); VOLUME(1)=DXCELL(I)*DYCELL(J)*DZCELL(K) ALLOCATE(XYZCEN(IAXIS:KAXIS,1:NCELL)); XYZCEN(IAXIS:KAXIS,1) = (/ XCELL(I),YCELL(J),ZCELL(K) /) - ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = .FALSE. + ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = NOT_BLOCKED ! Add one by one regular and gas cut faces: CT = 1; CCELEM(1,1) = 0 @@ -5266,7 +5281,7 @@ SUBROUTINE GET_REMAINING_CUTCELLS(NM) ALLOCATE(FACE_LIST(1:CC_NPARAM_CCFACE,1:NFACE_CELL)); FACE_LIST = CC_UNDEFINED ALLOCATE(VOLUME(1:NCELL)); VOLUME(1)=DXCELL(I)*DYCELL(J)*DZCELL(K) ALLOCATE(XYZCEN(IAXIS:KAXIS,1:NCELL)); XYZCEN(IAXIS:KAXIS,1) = (/ XCELL(I),YCELL(J),ZCELL(K) /) - ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = .FALSE. + ALLOCATE(NOADVANCE(1:NCELL)); NOADVANCE(1) = NOT_BLOCKED ! Add one by one regular and gas cut faces: CT = 1; CCELEM(1,1) = 0 @@ -6180,7 +6195,7 @@ SUBROUTINE BLOCK_CUT_CELL(NM,ICC,JCC,BLOCK_PHASE) LOWJ=M%CUT_FACE(IFCX)%CELL_LIST(2,HIGH_IND,JFCX) HIGJ=M%CUT_FACE(IFCX)%CELL_LIST(3,HIGH_IND,JFCX) IF(LOWI>0 .AND. LOWJ>0) THEN - IF(M%CUT_CELL(LOWI)%NOADVANCE(HIGI) .AND. M%CUT_CELL(LOWJ)%NOADVANCE(HIGJ)) CYCLE IFC_LOOP + IF(M%CUT_CELL(LOWI)%NOADVANCE(HIGI)>0 .AND. M%CUT_CELL(LOWJ)%NOADVANCE(HIGJ)>0) CYCLE IFC_LOOP ENDIF ENDIF @@ -6230,6 +6245,7 @@ SUBROUTINE BLOCK_CUT_CELL(NM,ICC,JCC,BLOCK_PHASE) ! Provide GEOM surface information to newly created INBOUNDARY face: M%CUT_FACE(IFC1)%BODTRI(1:2,JFC1) = (/ IBOD, ITRI /) M%CUT_FACE(IFC1)%SURF_INDEX(JFC1) = 0 ! Default surf. + M%CUT_FACE(IFC1)%CFACE_ORIGIN(JFC1) = M%CUT_CELL(ICC)%NOADVANCE(JCC) IF(IBOD>0) M%CUT_FACE(IFC1)%SURF_INDEX(JFC1) = GEOMETRY(IBOD)%SURFS(ITRI) M%CUT_FACE(IFC1)%NFACE = JFC1 ENDIF FACE_TYPE_IF @@ -6844,8 +6860,8 @@ SUBROUTINE BLOCK_CUT_CELL(NM,ICC,JCC,BLOCK_PHASE) LOWJ=M%CUT_FACE(IFCX)%CELL_LIST(2,HIGH_IND,JFCX) HIGJ=M%CUT_FACE(IFCX)%CELL_LIST(3,HIGH_IND,JFCX) IF(LOWI>0 .AND. LOWJ>0) THEN - IF(M%CUT_CELL(LOWI)%NOADVANCE(HIGI) .AND. & ! This is to drop this cut-face on the second hit. - M%CUT_CELL(LOWJ)%NOADVANCE(HIGJ) .AND. M%CUT_FACE(IFCX)%SHARED(JFCX)) THEN + IF(M%CUT_CELL(LOWI)%NOADVANCE(HIGI)>0 .AND. & ! This is to drop this cut-face on the second hit. + M%CUT_CELL(LOWJ)%NOADVANCE(HIGJ)>0 .AND. M%CUT_FACE(IFCX)%SHARED(JFCX)) THEN M%CUT_FACE(IFCX)%SHARED(JFCX) =.FALSE. CYCLE IFC_LOOP_2 ENDIF @@ -6874,7 +6890,7 @@ SUBROUTINE BLOCK_CUT_CELL(NM,ICC,JCC,BLOCK_PHASE) LOWJ=M%CUT_FACE(IFCX)%CELL_LIST(2,HIGH_IND,JFCX) HIGJ=M%CUT_FACE(IFCX)%CELL_LIST(3,HIGH_IND,JFCX) IF(LOWI>0 .AND. LOWJ>0) THEN - IF(M%CUT_CELL(LOWI)%NOADVANCE(HIGI) .AND. M%CUT_CELL(LOWJ)%NOADVANCE(HIGJ)) THEN + IF(M%CUT_CELL(LOWI)%NOADVANCE(HIGI)>0 .AND. M%CUT_CELL(LOWJ)%NOADVANCE(HIGJ)>0) THEN DROP_FACE=.TRUE. M%CUT_FACE(IFCX)%SHARED(JFCX) =.TRUE. ENDIF @@ -7858,7 +7874,7 @@ SUBROUTINE GET_CELL_LINK_INFO(NM) DO ICC=1,M%N_CUTCELL_MESH+M%N_GCCUTCELL_MESH M%CUT_CELL(ICC)%UNKZ(:) = CC_UNDEFINED DO JCC=1,M%CUT_CELL(ICC)%NCELL - IF (M%CUT_CELL(ICC)%NOADVANCE(JCC)) M%CUT_CELL(ICC)%IJK_LINK(1,JCC) = CC_SOLID + IF (M%CUT_CELL(ICC)%NOADVANCE(JCC)>0) M%CUT_CELL(ICC)%IJK_LINK(1,JCC) = CC_SOLID ENDDO ENDDO @@ -7873,7 +7889,7 @@ SUBROUTINE GET_CELL_LINK_INFO(NM) ! or cutcells inside an OBST. CCVOL_THRES = CCVOL_LINK * (M%DX(I)*M%DY(J)*M%DZ(K)) DO JCC=1,M%CUT_CELL(ICC)%NCELL - IF ( M%CUT_CELL(ICC)%NOADVANCE(JCC) ) CYCLE + IF ( M%CUT_CELL(ICC)%NOADVANCE(JCC)>0 ) CYCLE IF ( M%CUT_CELL(ICC)%VOLUME(JCC) > CCVOL_THRES) M%CUT_CELL(ICC)%UNKZ(JCC) = 1 ENDDO ENDIF @@ -7913,7 +7929,7 @@ SUBROUTINE GET_CELL_LINK_INFO(NM) CCVOL_THRES = CCVOL_LINK * (M%DX(I)*M%DY(J)*M%DZ(K)) JCC_LOOP_1 : DO JCC=1,CC%NCELL - IF ( CC%NOADVANCE(JCC) .OR. CC%UNKZ(JCC)>0 ) CYCLE + IF ( CC%NOADVANCE(JCC)>0 .OR. CC%UNKZ(JCC)>0 ) CYCLE CRTCELL_FLG = .FALSE. VAL_UNKZ = CC_UNDEFINED VAL_CVOL = CCVOL_THRES @@ -8019,7 +8035,7 @@ SUBROUTINE GET_CELL_LINK_INFO(NM) CCVOL_THRES = CCVOL_LINK * (M%DX(I)*M%DY(J)*M%DZ(K)) JCC_LOOP_2 : DO JCC=1,CC%NCELL - IF ( CC%NOADVANCE(JCC) .OR. CC%UNKZ(JCC)>0 ) CYCLE + IF ( CC%NOADVANCE(JCC)>0 .OR. CC%UNKZ(JCC)>0 ) CYCLE VAL_UNKZ = CC_UNDEFINED VAL_CVOL = -GEOMEPS @@ -8203,7 +8219,7 @@ SUBROUTINE GET_CELL_LINK_INFO(NM) ! JCC_LNK = JCC ! IF (JCC_LNK <= CC%NCELL) THEN ! DO JCC=1,CC%NCELL - ! IF ( CC%NOADVANCE(JCC) .OR. JCC==JCC_LNK ) CYCLE + ! IF ( CC%NOADVANCE(JCC)>0 .OR. JCC==JCC_LNK ) CYCLE ! CC%UNKZ(JCC) = CC%UNKZ(JCC_LNK) ! CC%IJK_LINK(1:KAXIS+2,JCC) = (/ CC_CUTCFE, I, J, K, JCC_LNK /) ! CC%LINK_LEV(JCC) = CC%LINK_LEV(JCC_LNK) - 1 @@ -8224,7 +8240,7 @@ SUBROUTINE GET_CELL_LINK_INFO(NM) CC => M%CUT_CELL(ICC); I = CC%IJK(IAXIS); J = CC%IJK(JAXIS); K = CC%IJK(KAXIS) IF (M%CELL(M%CELL_INDEX(I,J,K))%SOLID) CYCLE DO JCC=1,CC%NCELL - IF ( CC%NOADVANCE(JCC) .OR. CC%UNKZ(JCC)>0 ) CYCLE + IF ( CC%NOADVANCE(JCC)>0 .OR. CC%UNKZ(JCC)>0 ) CYCLE ULINK_COUNT = ULINK_COUNT + 1 ENDDO ENDDO @@ -8252,7 +8268,7 @@ SUBROUTINE GET_CELL_LINK_INFO(NM) CC => M%CUT_CELL(ICC); I = CC%IJK(IAXIS); J = CC%IJK(JAXIS); K = CC%IJK(KAXIS) IF (M%CELL(M%CELL_INDEX(I,J,K))%SOLID) CYCLE DO JCC=1,CC%NCELL - IF (CC%NOADVANCE(JCC) .OR. CC%UNKZ(JCC)>0) CYCLE + IF (CC%NOADVANCE(JCC)>0 .OR. CC%UNKZ(JCC)>0) CYCLE ULINK_COUNT = ULINK_COUNT + 1 WRITE(LU_UNLNK,'(I8,A,5I8,A,5F22.8)') & ULINK_COUNT,', I,J,K,ICC,JCC=',I,J,K,ICC,JCC,', X,Y,Z,CCVOL,CCVOL_CRT=',M%X(I),M%Y(J),M%Z(K), & @@ -8314,7 +8330,7 @@ SUBROUTINE GET_CELL_LINK_INFO(NM) DO ICC=1,M%N_CUTCELL_MESH DO JCC=1,M%CUT_CELL(ICC)%NCELL IF(M%CUT_CELL(ICC)%IJK_LINK(1,JCC)==CC_SOLID) THEN - M%CUT_CELL(ICC)%NOADVANCE(JCC) = .TRUE. + M%CUT_CELL(ICC)%NOADVANCE(JCC) = BLOCKED_UNLINK_CELL M%CUT_CELL(ICC)%LINK_LEV(JCC) = CC_UNDEFINED M%CUT_CELL(ICC)%IJK_LINK(2:5,JCC)= CC_UNDEFINED ELSEIF(M%CUT_CELL(ICC)%LINK_LEV(JCC) < M%FINEST_LINK_LEV) THEN @@ -17075,6 +17091,10 @@ SUBROUTINE FACE_REALLOC(NM,ICF,NVERT,NFACE,NSVERT,NSFACE,NVERTFACE_NEW) ALLOCATE(REAL1D(1:NFACE+NSFACE)); REAL1D = 1._EB CALL MOVE_ALLOC(FROM=REAL1D,TO=MESHES(NM)%CUT_FACE(ICF)%AREA_ADJUST) + ALLOCATE(INT1D(1:NFACE+NSFACE)); INT1D=NOT_BLOCKED + INT1D(1:NFACE) = MESHES(NM)%CUT_FACE(ICF)%CFACE_ORIGIN(1:NFACE) + CALL MOVE_ALLOC(FROM=INT1D,TO=MESHES(NM)%CUT_FACE(ICF)%CFACE_ORIGIN) + ENDIF ENDIF @@ -17131,6 +17151,7 @@ SUBROUTINE CUT_FACE_MOVE(CUT_FACE_FROM,CUT_FACE_TO) CALL MOVE_ALLOC(FROM=CUT_FACE_FROM%XYZCEN, TO=CUT_FACE_TO%XYZCEN) CALL MOVE_ALLOC(FROM=CUT_FACE_FROM%SHARED, TO=CUT_FACE_TO%SHARED) CALL MOVE_ALLOC(FROM=CUT_FACE_FROM%BLK_TAG, TO=CUT_FACE_TO%BLK_TAG) +CALL MOVE_ALLOC(FROM=CUT_FACE_FROM%CFACE_ORIGIN, TO=CUT_FACE_TO%CFACE_ORIGIN) CALL MOVE_ALLOC(FROM=CUT_FACE_FROM%LINK_LEV, TO=CUT_FACE_TO%LINK_LEV) CALL MOVE_ALLOC(FROM=CUT_FACE_FROM%INXAREA, TO=CUT_FACE_TO%INXAREA) CALL MOVE_ALLOC(FROM=CUT_FACE_FROM%INXSQAREA, TO=CUT_FACE_TO%INXSQAREA) @@ -17199,6 +17220,7 @@ SUBROUTINE FACE_DEALLOC(NM,ICF,DO_BNCF) ENDIF IF(ALLOCATED(MESHES(NM)%CUT_FACE(ICF)%SHARED)) DEALLOCATE(MESHES(NM)%CUT_FACE(ICF)%SHARED) IF(ALLOCATED(MESHES(NM)%CUT_FACE(ICF)%BLK_TAG)) DEALLOCATE(MESHES(NM)%CUT_FACE(ICF)%BLK_TAG) +IF(ALLOCATED(MESHES(NM)%CUT_FACE(ICF)%CFACE_ORIGIN)) DEALLOCATE(MESHES(NM)%CUT_FACE(ICF)%CFACE_ORIGIN) IF(ALLOCATED(MESHES(NM)%CUT_FACE(ICF)%LINK_LEV)) DEALLOCATE(MESHES(NM)%CUT_FACE(ICF)%LINK_LEV) IF(ALLOCATED(MESHES(NM)%CUT_FACE(ICF)%INXAREA)) DEALLOCATE(MESHES(NM)%CUT_FACE(ICF)%INXAREA) IF(ALLOCATED(MESHES(NM)%CUT_FACE(ICF)%INXSQAREA)) DEALLOCATE(MESHES(NM)%CUT_FACE(ICF)%INXSQAREA) @@ -17256,6 +17278,7 @@ SUBROUTINE NEW_FACE_ALLOC(NM,ICF,NVERT,NFACE,NVERTFACE,IBNDINT) IF(MESHES(NM)%CUT_FACE(ICF)%STATUS==CC_INBOUNDARY) THEN ALLOCATE(MESHES(NM)%CUT_FACE(ICF)%AREA_ADJUST(1:NFACE)); MESHES(NM)%CUT_FACE(ICF)%AREA_ADJUST = 1._EB + ALLOCATE(MESHES(NM)%CUT_FACE(ICF)%CFACE_ORIGIN(1:NFACE)); MESHES(NM)%CUT_FACE(ICF)%CFACE_ORIGIN = NOT_BLOCKED ELSE IF(PRESENT(IBNDINT)) THEN IF(IBNDINT>2) RETURN ! Gas cut-face not in block boundary. @@ -19905,7 +19928,7 @@ SUBROUTINE GET_CARTCELL_CUTCELLS(NM) INTEGER, ALLOCATABLE, DIMENSION(:,:) :: CCELEM REAL(EB),ALLOCATABLE, DIMENSION(:,:) :: XYZCEN REAL(EB),ALLOCATABLE, DIMENSION(:) :: VOL ! Cut-cell volumes. -LOGICAL, ALLOCATABLE, DIMENSION(:) :: NOADVANCE +INTEGER, ALLOCATABLE, DIMENSION(:) :: NOADVANCE INTEGER :: IFACE, IEDGE, ISEG, SEG(NOD1:NOD2), ICELL, NFACEI LOGICAL :: INLIST, TEST1, TEST2, NEWFACE @@ -20530,7 +20553,7 @@ SUBROUTINE GET_CARTCELL_CUTCELLS(NM) ENDIF CCELEM(1:NFACE_CELL+1,NCELL) = (/ NFACE_CELL, (IFACE, IFACE=1,NFACE_CELL) /) VOL(NCELL) = DXCELL(I)*DYCELL(J)*DZCELL(K) - NOADVANCE(NCELL) = .FALSE. + NOADVANCE(NCELL) = NOT_BLOCKED XYZCEN(IAXIS:KAXIS,NCELL) = (/ XCELL(I), YCELL(J), ZCELL(K) /) ELSE CYCLE_CELL_COND @@ -20565,7 +20588,7 @@ SUBROUTINE GET_CARTCELL_CUTCELLS(NM) ! Compute volumes and centroids for the found cut-cells: VOL(1:NCELL) = 0._EB - NOADVANCE(1:NCELL) = .FALSE. + NOADVANCE(1:NCELL) = NOT_BLOCKED XYZCEN(IAXIS:KAXIS,1:NCELL) = 0._EB DO ICELL=1,NCELL NP = CCELEM(1,ICELL) @@ -20578,8 +20601,9 @@ SUBROUTINE GET_CARTCELL_CUTCELLS(NM) ENDDO VOL(ICELL) = ABS(VOL(ICELL)) - ! Define if cut-cell is very small -> NOADVANCE(ICELL)=T: - IF(DO_NOADVANCE .AND. VOL(ICELL)/(DXCELL(I)*DYCELL(J)*DZCELL(K)) NOADVANCE(ICELL)=BLOCKED_SMALL_CELL: + IF(DO_NOADVANCE .AND. VOL(ICELL)/(DXCELL(I)*DYCELL(J)*DZCELL(K))DXCELL(I)*DYCELL(J)*DZCELL(K)) VOL(ICELL) = DXCELL(I)*DYCELL(J)*DZCELL(K) IF(VOL(ICELL) < GEOMEPS) THEN ! Volume too small for correct calculation of XYZCEN-> take cartcell centroid. IF(.NOT.DO_NOADVANCE .AND. VOL(ICELL)