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
62 changes: 40 additions & 22 deletions Source/ccib.f90
Original file line number Diff line number Diff line change
Expand Up @@ -271,7 +271,7 @@ SUBROUTINE GET_H_GUARD_CUTCELL(IPZ,HP)
! Local Variables:
INTEGER :: IW, II, JJ, KK, IIG, JJG, KKG, ICC, ICC_GHOST, ICC_INT, NOM, IIO, JJO, KKO, IOR, ICC_EXT
INTEGER :: X1AXIS, ICF_EXT, ISHF(IAXIS:KAXIS)
REAL(EB) :: H_INT, D_INT, D_GHOST, D_EXT, X_FACE, SUM_FLUX, A_INT, A_EXT, H_EXT, H_GHOST
REAL(EB) :: H_INT, D_INT, D_GHOST, D_EXT, X_FACE, XFN, SUM_FLUX, A_INT, A_EXT, H_EXT, H_GHOST
TYPE (WALL_TYPE), POINTER :: WC
TYPE (BOUNDARY_COORD_TYPE), POINTER :: BC
TYPE (EXTERNAL_WALL_TYPE), POINTER :: EWC
Expand Down Expand Up @@ -368,24 +368,30 @@ SUBROUTINE GET_H_GUARD_CUTCELL(IPZ,HP)
IF (PREDICTOR) THEN; H_EXT = OM%H(IIO,JJO,KKO)
ELSE; H_EXT = OM%HS(IIO,JJO,KKO)
ENDIF
! Neighbor seam-face index in NOM's own frame (opposite side of cell IIO,JJO,KKO).
ISHF(IAXIS:KAXIS) = 0; IF(IOR < 0) ISHF(X1AXIS) = -1
SELECT CASE(X1AXIS)
CASE(1)
A_EXT = MESHES(NOM)%DY(JJO) * MESHES(NOM)%DZ(KKO)
IF (ICC_EXT > 0) THEN; D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(IAXIS,1) - X_FACE)
ELSE; D_EXT = ABS(MESHES(NOM)%XC(IIO) - X_FACE)
A_EXT = MESHES(NOM)%DY(JJO) * MESHES(NOM)%DZ(KKO); XFN = MESHES(NOM)%X(IIO+ISHF(IAXIS))
IF (ICC_EXT > 0) THEN; D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(IAXIS,1) - XFN)
ELSE; D_EXT = ABS(MESHES(NOM)%XC(IIO) - XFN)
ENDIF
CASE(2)
A_EXT = MESHES(NOM)%DX(IIO) * MESHES(NOM)%DZ(KKO)
IF (ICC_EXT > 0) THEN; D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(JAXIS,1) - X_FACE)
ELSE; D_EXT = ABS(MESHES(NOM)%YC(JJO) - X_FACE)
A_EXT = MESHES(NOM)%DX(IIO) * MESHES(NOM)%DZ(KKO); XFN = MESHES(NOM)%Y(JJO+ISHF(JAXIS))
IF (ICC_EXT > 0) THEN; D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(JAXIS,1) - XFN)
ELSE; D_EXT = ABS(MESHES(NOM)%YC(JJO) - XFN)
ENDIF
CASE(3)
A_EXT = MESHES(NOM)%DX(IIO) * MESHES(NOM)%DY(JJO)
IF (ICC_EXT > 0) THEN; D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(KAXIS,1) - X_FACE)
ELSE; D_EXT = ABS(MESHES(NOM)%ZC(KKO) - X_FACE)
A_EXT = MESHES(NOM)%DX(IIO) * MESHES(NOM)%DY(JJO); XFN = MESHES(NOM)%Z(KKO+ISHF(KAXIS))
IF (ICC_EXT > 0) THEN; D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(KAXIS,1) - XFN)
ELSE; D_EXT = ABS(MESHES(NOM)%ZC(KKO) - XFN)
ENDIF
END SELECT
ISHF(IAXIS:KAXIS) = 0; IF(IOR < 0) ISHF(X1AXIS) = -1
! Refinement seam only: match GET_H_MATRIX_CC sub-face selection (skip solid/no-unknown partners).
IF (WC%BOUNDARY_TYPE == INTERPOLATED_BOUNDARY) THEN
IF (ALLOCATED(OM%MUNKH)) THEN; IF (OM%MUNKH(IIO,JJO,KKO) < 1) CYCLE; ENDIF
IF (MESHES(NOM)%FCVAR(IIO+ISHF(IAXIS),JJO+ISHF(JAXIS),KKO+ISHF(KAXIS),CC_FGSC,X1AXIS) == CC_SOLID) CYCLE
ENDIF
ICF_EXT = MESHES(NOM)%FCVAR(IIO+ISHF(IAXIS),JJO+ISHF(JAXIS),KKO+ISHF(KAXIS),CC_IDCF,X1AXIS)
IF (ICF_EXT > 0) A_EXT = SUM(MESHES(NOM)%CUT_FACE(ICF_EXT)%AREA(1:MESHES(NOM)%CUT_FACE(ICF_EXT)%NFACE))
IF (EWC%AREA_RATIO < 0.99_EB) A_EXT = A_INT
Expand Down Expand Up @@ -20067,7 +20073,7 @@ SUBROUTINE GET_H_MATRIX_CC(NM,NM1,IPZ)
INTEGER :: X1AXIS,IFACE,ICF,I,J,K,IND(LOW_IND:HIGH_IND),IND_LOC(LOW_IND:HIGH_IND),ILOC,JLOC,JCOL,IROW, &
LOCROW,LOCROW_1,LOCROW_2,IW,NOM,IIO,JJO,KKO,ICC_EXT,IND_INT,IND_EXT,ICFO,ISHF(IAXIS:KAXIS), &
ICC,II,JJ,KK,IIG,JJG,KKG
REAL(EB) :: AF,IDX,BIJ,KFACE(2,2),X_FACE,D_INT,D_EXT,AF_EXT
REAL(EB) :: AF,IDX,BIJ,KFACE(2,2),X_FACE,XFN,D_INT,D_EXT,AF_EXT
TYPE(CC_CUTFACE_TYPE), POINTER :: CF
TYPE(CC_RCFACE_TYPE), POINTER :: RCF
TYPE(EXTERNAL_WALL_TYPE), POINTER :: EWC
Expand Down Expand Up @@ -20187,24 +20193,30 @@ SUBROUTINE GET_H_MATRIX_CC(NM,NM1,IPZ)
DO IIO = EWC%IIO_MIN, EWC%IIO_MAX
ICC_EXT = MESHES(NOM)%CCVAR(IIO,JJO,KKO,CC_IDCC)
IND_EXT = OMESH(NOM)%MUNKH(IIO,JJO,KKO)
! Neighbor seam-face index in NOM's own frame (opposite side of cell IIO,JJO,KKO).
ISHF(IAXIS:KAXIS) = 0; IF(BOUNDARY_COORD(WC%BC_INDEX)%IOR < 0) ISHF(X1AXIS) = -1
SELECT CASE(X1AXIS)
CASE(IAXIS)
AF_EXT = MESHES(NOM)%DY(JJO) * MESHES(NOM)%DZ(KKO)
D_EXT = ABS(MESHES(NOM)%XC(IIO) - X_FACE)
XFN = MESHES(NOM)%X(IIO+ISHF(IAXIS)); D_EXT = ABS(MESHES(NOM)%XC(IIO) - XFN)
CASE(JAXIS)
AF_EXT = MESHES(NOM)%DX(IIO) * MESHES(NOM)%DZ(KKO)
D_EXT = ABS(MESHES(NOM)%YC(JJO) - X_FACE)
XFN = MESHES(NOM)%Y(JJO+ISHF(JAXIS)); D_EXT = ABS(MESHES(NOM)%YC(JJO) - XFN)
CASE(KAXIS)
AF_EXT = MESHES(NOM)%DX(IIO) * MESHES(NOM)%DY(JJO)
D_EXT = ABS(MESHES(NOM)%ZC(KKO) - X_FACE)
XFN = MESHES(NOM)%Z(KKO+ISHF(KAXIS)); D_EXT = ABS(MESHES(NOM)%ZC(KKO) - XFN)
END SELECT
IF (ICC_EXT > 0) THEN
! Note: IND_EXT comes from OMESH(NOM)%MUNKH which was populated
! via COPY_CC_MUNKH_TO_UNKH from communicated HS data
D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(X1AXIS,1) - X_FACE)
D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(X1AXIS,1) - XFN)
ENDIF
! Check for cut-face in NOM
ISHF(IAXIS:KAXIS) = 0; IF(BOUNDARY_COORD(WC%BC_INDEX)%IOR < 0) ISHF(X1AXIS) = -1
! Refinement seam only (see GET_H_GUARD_CUTCELL for PERIODIC note):
IF (WC%BOUNDARY_TYPE == INTERPOLATED_BOUNDARY) THEN
IF (IND_EXT < 1) CYCLE
IF (MESHES(NOM)%FCVAR(IIO+ISHF(IAXIS),JJO+ISHF(JAXIS),KKO+ISHF(KAXIS),CC_FGSC,X1AXIS) == CC_SOLID) CYCLE
ENDIF
ICFO = MESHES(NOM)%FCVAR(IIO+ISHF(IAXIS),JJO+ISHF(JAXIS),KKO+ISHF(KAXIS),CC_IDCF,X1AXIS)
IF (ICFO > 0) AF_EXT = SUM(MESHES(NOM)%CUT_FACE(ICFO)%AREA(1:MESHES(NOM)%CUT_FACE(ICFO)%NFACE))
IF (EWC%AREA_RATIO < 0.99_EB) AF_EXT = AF
Expand Down Expand Up @@ -20324,24 +20336,30 @@ SUBROUTINE GET_H_MATRIX_CC(NM,NM1,IPZ)
DO IIO = EWC%IIO_MIN, EWC%IIO_MAX
ICC_EXT = MESHES(NOM)%CCVAR(IIO,JJO,KKO,CC_IDCC)
IND_EXT = OMESH(NOM)%MUNKH(IIO,JJO,KKO)
! Neighbor seam-face index in NOM's own frame (opposite side of cell IIO,JJO,KKO).
ISHF(IAXIS:KAXIS) = 0; IF (BC%IOR < 0) ISHF(X1AXIS) = -1
SELECT CASE(X1AXIS)
CASE(IAXIS)
AF_EXT = MESHES(NOM)%DY(JJO) * MESHES(NOM)%DZ(KKO)
D_EXT = ABS(MESHES(NOM)%XC(IIO) - X_FACE)
XFN = MESHES(NOM)%X(IIO+ISHF(IAXIS)); D_EXT = ABS(MESHES(NOM)%XC(IIO) - XFN)
CASE(JAXIS)
AF_EXT = MESHES(NOM)%DX(IIO) * MESHES(NOM)%DZ(KKO)
D_EXT = ABS(MESHES(NOM)%YC(JJO) - X_FACE)
XFN = MESHES(NOM)%Y(JJO+ISHF(JAXIS)); D_EXT = ABS(MESHES(NOM)%YC(JJO) - XFN)
CASE(KAXIS)
AF_EXT = MESHES(NOM)%DX(IIO) * MESHES(NOM)%DY(JJO)
D_EXT = ABS(MESHES(NOM)%ZC(KKO) - X_FACE)
XFN = MESHES(NOM)%Z(KKO+ISHF(KAXIS)); D_EXT = ABS(MESHES(NOM)%ZC(KKO) - XFN)
END SELECT
IF (ICC_EXT > 0) THEN
! Note: IND_EXT comes from OMESH(NOM)%MUNKH which was populated
! via COPY_CC_MUNKH_TO_UNKH from communicated HS data
D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(X1AXIS,1) - X_FACE)
D_EXT = ABS(MESHES(NOM)%CUT_CELL(ICC_EXT)%XYZCEN(X1AXIS,1) - XFN)
ENDIF
! Check for cut-face in NOM
ISHF(IAXIS:KAXIS) = 0; IF (BC%IOR < 0) ISHF(X1AXIS) = -1
! Refinement seam only (see GET_H_GUARD_CUTCELL for PERIODIC note):
IF (WC%BOUNDARY_TYPE == INTERPOLATED_BOUNDARY) THEN
IF (IND_EXT < 1) CYCLE
IF (MESHES(NOM)%FCVAR(IIO+ISHF(IAXIS),JJO+ISHF(JAXIS),KKO+ISHF(KAXIS),CC_FGSC,X1AXIS) == CC_SOLID) CYCLE
ENDIF
ICFO = MESHES(NOM)%FCVAR(IIO+ISHF(IAXIS),JJO+ISHF(JAXIS),KKO+ISHF(KAXIS),CC_IDCF,X1AXIS)
IF (ICFO > 0) AF_EXT = SUM(MESHES(NOM)%CUT_FACE(ICFO)%AREA(1:MESHES(NOM)%CUT_FACE(ICFO)%NFACE))
IF (EWC%AREA_RATIO < 0.99_EB) AF_EXT = AF
Expand Down
33 changes: 32 additions & 1 deletion Source/geom.f90
Original file line number Diff line number Diff line change
Expand Up @@ -6970,6 +6970,26 @@ SUBROUTINE CUT_CELL_FACE_ARRAYS_CLEANUP(NM)
ENDDO
ENDDO

! Invalidate CCVAR at (I,J,K) pointing at cut-cells dropped from the active pool:
DO K=-CCGUARD,M%KBAR+CCGUARD
DO J=-CCGUARD,M%JBAR+CCGUARD
DO I=-CCGUARD,M%IBAR+CCGUARD
ICC = M%CCVAR(I,J,K,CC_IDCC)
IF (ICC < 1) CYCLE
IF (ICC > M%N_CUTCELL_MESH+M%N_GCCUTCELL_MESH) THEN
M%CCVAR(I,J,K,CC_CGSC) = CC_SOLID
M%CCVAR(I,J,K,CC_IDCC) = CC_UNDEFINED
ELSEIF (M%CUT_CELL(ICC)%NCELL < 1 .OR. &
M%CUT_CELL(ICC)%IJK(IAXIS) /= I .OR. &
M%CUT_CELL(ICC)%IJK(JAXIS) /= J .OR. &
M%CUT_CELL(ICC)%IJK(KAXIS) /= K) THEN
M%CCVAR(I,J,K,CC_CGSC) = CC_SOLID
M%CCVAR(I,J,K,CC_IDCC) = CC_UNDEFINED
ENDIF
ENDDO
ENDDO
ENDDO

DEALLOCATE(CCIND,CFIND)

RETURN
Expand Down Expand Up @@ -21935,9 +21955,20 @@ END SUBROUTINE CUT_CELL_MOVE
SUBROUTINE CELL_DEALLOC(NM,ICC)

INTEGER, INTENT(IN) :: NM,ICC
INTEGER :: I,J,K

MESHES(NM)%CUT_CELL(ICC)%NCELL = 0
IF (.NOT.ALLOCATED(MESHES(NM)%CUT_CELL(ICC)%CCELEM)) RETURN
IF (.NOT.ALLOCATED(MESHES(NM)%CUT_CELL(ICC)%CCELEM)) THEN
I = MESHES(NM)%CUT_CELL(ICC)%IJK(IAXIS)
J = MESHES(NM)%CUT_CELL(ICC)%IJK(JAXIS)
K = MESHES(NM)%CUT_CELL(ICC)%IJK(KAXIS)
IF (I>=LBOUND(MESHES(NM)%CCVAR,1) .AND. I<=UBOUND(MESHES(NM)%CCVAR,1) .AND. &
J>=LBOUND(MESHES(NM)%CCVAR,2) .AND. J<=UBOUND(MESHES(NM)%CCVAR,2) .AND. &
K>=LBOUND(MESHES(NM)%CCVAR,3) .AND. K<=UBOUND(MESHES(NM)%CCVAR,3)) THEN
IF (MESHES(NM)%CCVAR(I,J,K,CC_IDCC)==ICC) MESHES(NM)%CCVAR(I,J,K,CC_IDCC) = CC_UNDEFINED
ENDIF
RETURN
ENDIF

! Deallocate ICC entries:
DEALLOCATE(MESHES(NM)%CUT_CELL(ICC)%CCELEM)
Expand Down
Loading
Loading