diff --git a/Manuals/FDS_User_Guide/FDS_User_Guide.tex b/Manuals/FDS_User_Guide/FDS_User_Guide.tex index 22f3fc34537..44c1288be5b 100644 --- a/Manuals/FDS_User_Guide/FDS_User_Guide.tex +++ b/Manuals/FDS_User_Guide/FDS_User_Guide.tex @@ -2445,7 +2445,7 @@ \subsubsection{Limitations} \begin{enumerate} \item {\ct HT3D} cannot be applied to an exterior boundary; it must be applied to an {\ct OBST} for which at least one face in each coordinate direction is exposed. \item Avoid contact between 3-D and 1-D solids. If two sides of a 3-D solid touch 1-D solids, there will be no lateral heat conduction computed in that particular direction. -\item If your 3-D obstruction extends beyond meshes that abut, add the parameter {\ct NEIGHBOR\_SEPARATION\_DISTANCE} to the {\ct MISC} line. Any mesh within this distance of another mesh will share geometry information for use in the 3-D heat conduction calculation. +\item If your 3-D obstruction extends beyond meshes that abut, FDS uses a special algorithm to identify all the meshes where this obstruction lives, and also those obstructions connected to it. However, if this algorithm fails to detect all the meshes and you recieve an error stating that there is a problem, add the parameter {\ct NEIGHBOR\_SEPARATION\_DISTANCE} to the {\ct MISC} line. Any mesh within this distance of another mesh will share geometry information for use in the 3-D heat conduction calculation. If you set this parameter to a value larger than the width of the computational domain, then all meshes will establish communication channels for exchanging boundary information. \item By default, the interior nodes are clustered near the surface and stretched out deeper within the solid. If you want to maintain uniform spacing, set {\ct CELL\_SIZE} on the {\ct SURF} line to indicate the desired interior node spacing. The {\ct CELL\_SIZE} is typically the same as the gas phase cells. \end{enumerate} diff --git a/Source/read.f90 b/Source/read.f90 index f9a12dedbd1..425fdb7b13d 100644 --- a/Source/read.f90 +++ b/Source/read.f90 @@ -49,6 +49,14 @@ MODULE READ_INPUT TYPE(MATERIAL_TYPE), POINTER :: ML=>NULL() TYPE(REACTION_TYPE), POINTER :: RN=>NULL() LOGICAL :: RETURN_BEFORE_STOP_FILE=.FALSE., RETURN_BEFORE_SIM_MODE=.FALSE. +INTEGER :: N_HT3D_SURF_LINES=0,N_HT3D_OBST=0 +CHARACTER(LABEL_LENGTH), DIMENSION(20) :: HT3D_SURF_LIST='null' +TYPE HT3D_OBST_TYPE + INTEGER :: GROUP_INDEX=0 + REAL(EB) :: XS,XF,YS,YF,ZS,ZF +END TYPE +TYPE(HT3D_OBST_TYPE), DIMENSION(:), ALLOCATABLE, TARGET :: HT3D_OBST + CONTAINS @@ -103,6 +111,7 @@ SUBROUTINE READ_DATA(DT) CALL READ_MOVE ; CALL CHECK_STOP_STATUS ; IF (STOP_STATUS/=NO_STOP) RETURN CALL READ_MULT ; CALL CHECK_STOP_STATUS ; IF (STOP_STATUS/=NO_STOP) RETURN CALL READ_MESH ; CALL CHECK_STOP_STATUS ; IF (STOP_STATUS/=NO_STOP) RETURN +CALL PROC_MESH ; CALL CHECK_STOP_STATUS ; IF (STOP_STATUS/=NO_STOP) RETURN CALL READ_TRAN ; CALL CHECK_STOP_STATUS ; IF (STOP_STATUS/=NO_STOP) RETURN CALL READ_TIME(DT); CALL CHECK_STOP_STATUS ; IF (STOP_STATUS/=NO_STOP) RETURN CALL READ_PRES ; CALL CHECK_STOP_STATUS ; IF (STOP_STATUS/=NO_STOP) RETURN @@ -539,15 +548,13 @@ END FUNCTION READ_STOP SUBROUTINE READ_MESH -USE COMP_FUNCTIONS, ONLY : SEARCH_INPUT_FILE -INTEGER :: IJK(3),NM,NM2,CURRENT_MPI_PROCESS,MPI_PROCESS,RGB(3),N_MESH_NEW,N,II,JJ,KK,NMESHES_READ,NNN,JBAR_OLD_VALUE -INTEGER, ALLOCATABLE, DIMENSION(:) :: NEIGHBOR_LIST -LOGICAL :: OVERLAPPING_X,OVERLAPPING_Y,OVERLAPPING_Z,POSSIBLY_PERIODIC,CYLINDRICAL_OLD_VALUE,PERIODIC_FOUND_IN_FILE +INTEGER :: IJK(3),NM,CURRENT_MPI_PROCESS,MPI_PROCESS,RGB(3),N_MESH_NEW,N,II,JJ,KK,NMESHES_READ,NNN,JBAR_OLD_VALUE +LOGICAL :: CYLINDRICAL_OLD_VALUE REAL(EB) :: XB1,XB2,XB3,XB4,XB5,XB6 CHARACTER(25) :: COLOR CHARACTER(LABEL_LENGTH) :: MULT_ID,TRNX_ID,TRNY_ID,TRNZ_ID NAMELIST /MESH/ CHECK_MESH_ALIGNMENT,COLOR,CYLINDRICAL,FYI,ID,IJK,MPI_PROCESS,MULT_ID,RGB,TRNX_ID,TRNY_ID,TRNZ_ID,XB -TYPE (MESH_TYPE), POINTER :: M,M2 +TYPE (MESH_TYPE), POINTER :: M TYPE (MULTIPLIER_TYPE), POINTER :: MR NMESHES = 0 @@ -842,29 +849,117 @@ SUBROUTINE READ_MESH ENDIF ENDDO -! Determine mesh neighbors. MESH_SEPARATION_DISTANCE is a very small length -! used to determine if there are periodic boundaries. NEIGHBOR_SEPARATION_DISANCE -! is the distance beyond which no information or message passing is assumed -! between the meshes. Its value is deliberately complicated to avoid having two -! meshes separated by exactly that distance. +END SUBROUTINE READ_MESH + + +!> \brief Determine mesh neighborhoods for MPI communications + +SUBROUTINE PROC_MESH + +USE COMP_FUNCTIONS, ONLY : SEARCH_INPUT_FILE +INTEGER :: NM,NM2,I,II,III,N_GROUPS=0 +INTEGER, ALLOCATABLE, DIMENSION(:) :: NEIGHBOR_LIST +LOGICAL :: OVERLAPPING_X,OVERLAPPING_Y,OVERLAPPING_Z,POSSIBLY_PERIODIC,PERIODIC_FOUND_IN_FILE +TYPE (MESH_TYPE), POINTER :: M,M2 +TYPE (HT3D_OBST_TYPE), POINTER :: HO,HO2,HO3 +TYPE MESH_COMM_TYPE + INTEGER :: N_MESHES=0 + INTEGER, DIMENSION(100) :: LIST=0 +END TYPE MESH_COMM_TYPE +TYPE(MESH_COMM_TYPE), ALLOCATABLE, DIMENSION(:), TARGET :: MESH_COMM +TYPE(MESH_COMM_TYPE), POINTER :: MC + +! Read the SURF lines and make a list, HT3D_SURF_LIST, of those that are HT3D or VARIABLE_THICKNESS + +CALL READ_SURF(QUICK_READ=.TRUE.) + +! If there are HT3D solids, determine the indices of the meshes that contain connected HT3D OBSTs + +IF_HT3D: IF (N_HT3D_SURF_LINES>0) THEN + + ! Read the OBST lines and make a list, HT3D_OBST, of those that have HT3D or VARIABLE_THICKNESS SURF lines + + ALLOCATE(HT3D_OBST(100)) ! Initial allocation -- this array can be reallocated if needed + CALL READ_OBST(QUICK_READ=.TRUE.) + + ! Assign each HT3D_OBST an integer, GROUP_INDEX, which indicates which group of obstructions it is connected to. + ! N_GROUPS is the number of these blocks of connected OBSTs. + + N_GROUPS = 1 + + HT3D_OBST_LOOP: DO I=1,N_HT3D_OBST + HO => HT3D_OBST(I) + IF (I==1) THEN + HO%GROUP_INDEX = 1 + ELSE + HO%GROUP_INDEX = 1000000 + DO II=I-1,1,-1 + HO2 => HT3D_OBST(II) + IF (HO%XS>HO2%XF.OR.HO%XFHO2%YF.OR.HO%YFHO2%ZF.OR.HO%ZFHO%GROUP_INDEX) THEN ! If current OBST is connected to one already assigned a GROUP_INDEX, + ! and that GROUP_INDEX is greater than that of the OBST, reduce the + ! N_GROUPS and GROUP_INDEXs of OBSTs with higher values + N_GROUPS = N_GROUPS - 1 + DO III=1,I-1 + HO3 => HT3D_OBST(III) + IF (HO3%GROUP_INDEX>HO%GROUP_INDEX) HO3%GROUP_INDEX = HO3%GROUP_INDEX - 1 + ENDDO + ENDIF + HO%GROUP_INDEX = MIN(HO%GROUP_INDEX,HO2%GROUP_INDEX) + IF (HO%GROUP_INDEX==N_GROUPS) CYCLE HT3D_OBST_LOOP + ENDDO + IF (HO%GROUP_INDEX>100000) THEN ! The current OBST is not connected to any previous OBST. Assign it a new GROUP_INDEX. + N_GROUPS = N_GROUPS + 1 + HO%GROUP_INDEX = N_GROUPS + ENDIF + ENDIF + ENDDO HT3D_OBST_LOOP + + ! Now determine all the MESH indices corresponding to each group of connected OBSTs + + ALLOCATE(MESH_COMM(N_GROUPS)) + + DO I=1,N_HT3D_OBST + HO => HT3D_OBST(I) + MC => MESH_COMM(HO%GROUP_INDEX) + MESH_LOOP: DO NM=1,NMESHES + M => MESHES(NM) + IF (HO%XS>M%XF.OR.HO%XFM%YF.OR.HO%YFM%ZF.OR.HO%ZF MESHES(NM) M%N_NEIGHBORING_MESHES = 0 NEIGHBOR_LIST = 0 + + ! Add adjacent meshes to the neighborhood of MESH NM + DO NM2=1,NMESHES M2 => MESHES(NM2) OVERLAPPING_X = .TRUE. @@ -888,16 +983,31 @@ SUBROUTINE READ_MESH M%N_NEIGHBORING_MESHES = M%N_NEIGHBORING_MESHES + 1 NEIGHBOR_LIST(M%N_NEIGHBORING_MESHES) = NM2 ENDDO + + ! Add meshes containing the HT3D_OBST groups to the neighborhood of MESH NM + + DO I=1,N_GROUPS + MC => MESH_COMM(I) + IF (ANY(MC%LIST==NM,DIM=1)) THEN + DO II=1,MC%N_MESHES + IF (ANY(NEIGHBOR_LIST==MC%LIST(II))) CYCLE + M%N_NEIGHBORING_MESHES = M%N_NEIGHBORING_MESHES + 1 + NEIGHBOR_LIST(M%N_NEIGHBORING_MESHES) = MC%LIST(II) + ENDDO + ENDIF + ENDDO + + ! Save the list of neighboring meshes into an array + ALLOCATE(M%NEIGHBORING_MESH(M%N_NEIGHBORING_MESHES)) DO I=1,M%N_NEIGHBORING_MESHES M%NEIGHBORING_MESH(I) = NEIGHBOR_LIST(I) ENDDO ENDDO -DEALLOCATE(NEIGHBOR_LIST) -REWIND(LU_INPUT) ; INPUT_FILE_LINE_NUMBER = 0 +DEALLOCATE(NEIGHBOR_LIST) -END SUBROUTINE READ_MESH +END SUBROUTINE PROC_MESH !> \brief Read the TRAN namelist lines and compute the polynomial transform function for the vertical coordinate @@ -7051,10 +7161,11 @@ END SUBROUTINE PROC_MATL !> \brief Read the SURF namelist lines -SUBROUTINE READ_SURF +SUBROUTINE READ_SURF(QUICK_READ) USE MATH_FUNCTIONS, ONLY : GET_RAMP_INDEX USE DEVICE_VARIABLES, ONLY : PROPERTY_TYPE +LOGICAL, INTENT(IN), OPTIONAL :: QUICK_READ CHARACTER(LABEL_LENGTH) :: PART_ID,RAMP_MF(MAX_SPECIES),RAMP_Q,RAMP_V,RAMP_T,RAMP_T_I,MATL_ID(MAX_LAYERS,MAX_MATERIALS),& PROFILE,BACKING,GEOMETRY,RAMP_EF,RAMP_PART,NAME_LIST(MAX_MATERIALS*MAX_LAYERS),& SPEC_ID(MAX_SPECIES),RAMP_TMP_BACK,RAMP_TMP_GAS_BACK,RAMP_TMP_GAS_FRONT,& @@ -7130,11 +7241,17 @@ SUBROUTINE READ_SURF CALL CHECKREAD('SURF',LU_INPUT,IOS) ; IF (STOP_STATUS==SETUP_STOP) RETURN IF (IOS==1) EXIT COUNT_SURF_LOOP ID = 'null' + HT3D = .FALSE. + VARIABLE_THICKNESS = .FALSE. READ(LU_INPUT,SURF,ERR=34,IOSTAT=IOS) IF (ID=='null') THEN WRITE(MESSAGE,'(A)') 'ERROR: SURF line must have an ID' CALL SHUTDOWN(MESSAGE) ; RETURN ENDIF + IF (PRESENT(QUICK_READ) .AND. (HT3D .OR. VARIABLE_THICKNESS)) THEN + N_HT3D_SURF_LINES = N_HT3D_SURF_LINES + 1 + HT3D_SURF_LIST(N_HT3D_SURF_LINES) = ID + ENDIF N_SURF = N_SURF + 1 34 IF (IOS>0) THEN WRITE(MESSAGE,'(A,I0,A,I0)') 'ERROR: Problem with SURF number ',N_SURF+1,', line number ',INPUT_FILE_LINE_NUMBER @@ -7142,6 +7259,10 @@ SUBROUTINE READ_SURF ENDIF ENDDO COUNT_SURF_LOOP +! Special case where SURF lines are scanned, looking for presence of HT3D + +IF (PRESENT(QUICK_READ)) RETURN + ! Allocate the SURFACE derived type, leaving space for SURF entries not defined explicitly by the user N_SURF_RESERVED = 9 @@ -9532,16 +9653,17 @@ END SUBROUTINE READ_TABL !> \brief Read the OBSTruction namelist lines -SUBROUTINE READ_OBST +SUBROUTINE READ_OBST(QUICK_READ) USE MEMORY_FUNCTIONS, ONLY: REALLOCATE_CELL USE GEOMETRY_FUNCTIONS, ONLY: BLOCK_CELL,CIRCLE_CELL_INTERSECTION_AREA,SEARCH_OTHER_MESHES USE COMPLEX_GEOMETRY, ONLY: INTERSECT_CONE_AABB,INTERSECT_CYLINDER_AABB,INTERSECT_SPHERE_AABB,INTERSECT_OBB_AABB,ROTATION_MATRIX USE MISC_FUNCTIONS, ONLY: PROCESS_MESH_NEIGHBORHOOD +LOGICAL, INTENT(IN), OPTIONAL :: QUICK_READ TYPE(OBSTRUCTION_TYPE), POINTER :: OB2,OBT TYPE(MULTIPLIER_TYPE), POINTER :: MR TYPE(OBSTRUCTION_TYPE), DIMENSION(:), ALLOCATABLE, TARGET :: TEMP_OBSTRUCTION -INTEGER :: NM,NOM,NOM2,N_OBST_O,IC,N,NN,NNN,NNNN,N_NEW_OBST,RGB(3),N_OBST_DIM,II,JJ,KK,MULT_INDEX,SHAPE_TYPE,IIO,JJO,KKO,IOR +INTEGER :: I,NM,NOM,NOM2,N_OBST_O,IC,N,NN,NNN,NNNN,N_NEW_OBST,RGB(3),N_OBST_DIM,II,JJ,KK,MULT_INDEX,SHAPE_TYPE,IIO,JJO,KKO,IOR CHARACTER(LABEL_LENGTH) :: ID,DEVC_ID,SHAPE,SURF_ID,SURF_ID_INTERIOR,SURF_IDS(3),SURF_ID6(6),CTRL_ID,MULT_ID,MATL_ID(MAX_MATERIALS) CHARACTER(25) :: COLOR LOGICAL :: OVERLAY,IS_INTERSECT @@ -9556,11 +9678,11 @@ SUBROUTINE READ_OBST MESH_LOOP: DO NM=1,NMESHES - M => MESHES(NM) - - IF (.NOT.PROCESS_MESH_NEIGHBORHOOD(NM)) CYCLE MESH_LOOP - - CALL POINT_TO_MESH(NM) + IF (.NOT.PRESENT(QUICK_READ)) THEN + M => MESHES(NM) + IF (.NOT.PROCESS_MESH_NEIGHBORHOOD(NM)) CYCLE MESH_LOOP + CALL POINT_TO_MESH(NM) + ENDIF ! Count OBST lines @@ -9571,6 +9693,9 @@ SUBROUTINE READ_OBST CALL CHECKREAD('OBST',LU_INPUT,IOS) ; IF (STOP_STATUS==SETUP_STOP) RETURN IF (IOS==1) EXIT COUNT_OBST_LOOP MULT_ID = 'null' + SURF_ID = 'null' + SURF_IDS = 'null' + SURF_ID6 = 'null' READ(LU_INPUT,NML=OBST,END=1,ERR=2,IOSTAT=IOS) CALL CHECK_XB(XB) MULT_INDEX = -1 @@ -9609,6 +9734,22 @@ SUBROUTINE READ_OBST XB6 = XB(6) + MR%DZ0 + II*MR%DXB(6) ENDIF N_OBST_O = N_OBST_O + 1 + IF (PRESENT(QUICK_READ)) THEN + DO I=1,N_HT3D_SURF_LINES + IF (HT3D_SURF_LIST(I)==SURF_ID .OR. ANY(HT3D_SURF_LIST(I)==SURF_IDS) .OR. & + ANY(HT3D_SURF_LIST(I)==SURF_ID6)) THEN + IF (SIZE(HT3D_OBST)<=N_HT3D_OBST) CALL REALLOCATE_HT3D_OBST(SIZE(HT3D_OBST),SIZE(HT3D_OBST)+100) + N_HT3D_OBST = N_HT3D_OBST + 1 + HT3D_OBST(N_HT3D_OBST)%XS = XB1 + HT3D_OBST(N_HT3D_OBST)%XF = XB2 + HT3D_OBST(N_HT3D_OBST)%YS = XB3 + HT3D_OBST(N_HT3D_OBST)%YF = XB4 + HT3D_OBST(N_HT3D_OBST)%ZS = XB5 + HT3D_OBST(N_HT3D_OBST)%ZF = XB6 + ENDIF + ENDDO + CYCLE I_MULT_LOOP2 + ENDIF IF (XB1>M%XF+M%DX(M%IBAR) .OR. XB2M%YF+M%DY(M%JBAR) .OR. XB4M%ZF+M%DZ(M%KBAR) .OR. XB6