FDS Source: Issue firemodels#5866. Add routine to find HT3D wall cells.

Source/func.f90

@@ -1105,6 +1105,107 @@ SUBROUTINE BLOCK_CELL(NM,I1,I2,J1,J2,K1,K2,IVAL,OBST_INDEX)
 END SUBROUTINE BLOCK_CELL
 
 
+SUBROUTINE ASSIGN_HT3D_WALL_INDICES(NM)
+
+! For each cell with CELL_INDEX=IC in an HT3D solid obstruction, fill in the array MESHES(NM)%WALL_INDEX(IC,-3:3) such that:
+! WALL_INDEX(IC,-3) is the wall index of the bottom of the solid obstruction
+! WALL_INDEX(IC,-2) is the wall index of the back   of the solid obstruction
+! WALL_INDEX(IC,-1) is the wall index of the left   of the solid obstruction
+! WALL_INDEX(IC, 1) is the wall index of the right  of the solid obstruction
+! WALL_INDEX(IC, 2) is the wall index of the front  of the solid obstruction
+! WALL_INDEX(IC, 3) is the wall index of the top    of the solid obstruction
+! WALL_INDEX(IC, 0) is the wall index of the nearest surface of the solid obstruction
+
+INTEGER, INTENT(IN) :: NM
+INTEGER :: I,J,K,N,IC,II,JJ,KK,ICN,CELL_COUNT(-3:3)
+TYPE(MESH_TYPE), POINTER :: M
+TYPE(OBSTRUCTION_TYPE), POINTER :: OB
+
+M => MESHES(NM)
+
+DO N=1,M%N_OBST
+   OB => M%OBSTRUCTION(N)
+   IF (OB%HT3D) THEN
+      K_LOOP: DO K=OB%K1+1,OB%K2
+         J_LOOP: DO J=OB%J1+1,OB%J2
+            I_LOOP: DO I=OB%I1+1,OB%I2
+               IC = M%CELL_INDEX(I,J,K)
+               IF (.NOT.M%SOLID(IC)) CYCLE I_LOOP
+               M%WALL_INDEX(IC,:) = 0
+               CELL_COUNT = 0
+               CELL_COUNT(0) = 1000000
+
+               MARCH_RIGHT: DO II=I+1,M%IBAR
+                  ICN = M%CELL_INDEX(II,J,K)
+                  CELL_COUNT(1) = CELL_COUNT(1) + 1
+                  IF (.NOT.M%SOLID(ICN)) THEN
+                     M%WALL_INDEX(IC,1) = M%WALL_INDEX(ICN,-1)
+                     EXIT MARCH_RIGHT
+                  ENDIF
+                  IF (II==M%IBAR) CELL_COUNT(1) = 1000000
+               ENDDO MARCH_RIGHT
+
+               MARCH_LEFT: DO II=I-1,1,-1
+                  ICN = M%CELL_INDEX(II,J,K)
+                  CELL_COUNT(-1) = CELL_COUNT(-1) + 1
+                  IF (.NOT.M%SOLID(ICN)) THEN
+                     M%WALL_INDEX(IC,-1) = M%WALL_INDEX(ICN,1)
+                     EXIT MARCH_LEFT
+                  ENDIF
+                  IF (II==1) CELL_COUNT(-1) = 1000000
+               ENDDO MARCH_LEFT
+
+               MARCH_FORWARD: DO JJ=J+1,M%JBAR
+                  ICN = M%CELL_INDEX(I,JJ,K)
+                  CELL_COUNT(2) = CELL_COUNT(2) + 1
+                  IF (.NOT.M%SOLID(ICN)) THEN
+                     M%WALL_INDEX(IC,2) = M%WALL_INDEX(ICN,-2)
+                     EXIT MARCH_FORWARD
+                  ENDIF
+                  IF (JJ==M%JBAR) CELL_COUNT(2) = 1000000
+               ENDDO MARCH_FORWARD
+
+               MARCH_BACK: DO JJ=J-1,1,-1
+                  ICN = M%CELL_INDEX(I,JJ,K)
+                  CELL_COUNT(-2) = CELL_COUNT(-2) + 1
+                  IF (.NOT.M%SOLID(ICN)) THEN
+                     M%WALL_INDEX(IC,-2) = M%WALL_INDEX(ICN,2)
+                     EXIT MARCH_BACK
+                  ENDIF
+                  IF (JJ==1) CELL_COUNT(-2) = 1000000
+               ENDDO MARCH_BACK
+
+               MARCH_UP: DO KK=K+1,M%KBAR
+                  ICN = M%CELL_INDEX(I,J,KK)
+                  CELL_COUNT(3) = CELL_COUNT(3) + 1
+                  IF (.NOT.M%SOLID(ICN)) THEN
+                     M%WALL_INDEX(IC,3) = M%WALL_INDEX(ICN,-3)
+                     EXIT MARCH_UP
+                  ENDIF
+                  IF (KK==M%KBAR) CELL_COUNT(3) = 1000000
+               ENDDO MARCH_UP
+
+               MARCH_DOWN: DO KK=K-1,1,-1
+                  ICN = M%CELL_INDEX(I,J,KK)
+                  CELL_COUNT(-3) = CELL_COUNT(-3) + 1
+                  IF (.NOT.M%SOLID(ICN)) THEN
+                     M%WALL_INDEX(IC,-3) = M%WALL_INDEX(ICN,3)
+                     EXIT MARCH_DOWN
+                  ENDIF
+                  IF (KK==1) CELL_COUNT(-3) = 1000000
+               ENDDO MARCH_DOWN
+
+               M%WALL_INDEX(IC,0) = M%WALL_INDEX(IC,MINLOC(CELL_COUNT,DIM=1)-4)
+
+            ENDDO I_LOOP
+         ENDDO J_LOOP
+      ENDDO K_LOOP
+   ENDIF
+ENDDO
+
+END SUBROUTINE ASSIGN_HT3D_WALL_INDICES
+
+
 LOGICAL FUNCTION INTERIOR(XX,YY,ZZ)
 
 INTEGER NM

Source/init.f90

@@ -26,7 +26,7 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM)
 USE PHYSICAL_FUNCTIONS, ONLY: GET_VISCOSITY,GET_SPECIFIC_GAS_CONSTANT,GET_SPECIFIC_HEAT,LES_FILTER_WIDTH_FUNCTION
 USE MEMORY_FUNCTIONS, ONLY: COMPUTE_ONE_D_STORAGE_DIMENSIONS,COMPUTE_PARTICLE_STORAGE_DIMENSIONS,COMPUTE_WALL_STORAGE_DIMENSIONS,&
 COMPUTE_CFACE_STORAGE_DIMENSIONS
-USE GEOMETRY_FUNCTIONS, ONLY: ASSIGN_PRESSURE_ZONE
+USE GEOMETRY_FUNCTIONS, ONLY: ASSIGN_PRESSURE_ZONE,ASSIGN_HT3D_WALL_INDICES
 USE RADCONS, ONLY: UIIDIM
 USE CONTROL_VARIABLES
 USE MATH_FUNCTIONS, ONLY:EVALUATE_RAMP
@@ -984,6 +984,10 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM)
       IF (WC%BOUNDARY_TYPE==OPEN_BOUNDARY .AND. MEAN_FORCING(ABS(IOR))) M%APPLY_SPONGE_LAYER(IOR) = .TRUE.
    ENDDO
 ENDIF
+
+! Find exterior wall cells of HT3D solids
+
+CALL ASSIGN_HT3D_WALL_INDICES(NM)
 
 END SUBROUTINE INITIALIZE_MESH_VARIABLES_1