diff --git a/Source/dump.f90 b/Source/dump.f90 index e194dde29c..6c78f575a4 100644 --- a/Source/dump.f90 +++ b/Source/dump.f90 @@ -2223,12 +2223,13 @@ SUBROUTINE WRITE_SMOKEVIEW_FILE DO N=1,M%N_OBST OB => M%OBSTRUCTION(N) IF (OB%COLOR_INDICATOR/=-3) THEN - WRITE(MYSTR,'(8I5,A,L1)') OB%I1,OB%I2,OB%J1,OB%J2,OB%K1,OB%K2,OB%COLOR_INDICATOR,OB%TYPE_INDICATOR, & - ' ! ',OB%REMOVABLE; CALL ADDSTR + WRITE(MYSTR,'(8I5,A,L1,1X,6I2)') OB%I1,OB%I2,OB%J1,OB%J2,OB%K1,OB%K2,OB%COLOR_INDICATOR,OB%TYPE_INDICATOR, & + ' ! ',OB%REMOVABLE,OB%EXPOSED_FACE_INDEX(1:6) ELSE - WRITE(MYSTR,'(8I5,4F13.5,A,L1)') OB%I1,OB%I2,OB%J1,OB%J2,OB%K1,OB%K2,OB%COLOR_INDICATOR,OB%TYPE_INDICATOR, & - REAL(OB%RGB,FB)/255._FB, OB%TRANSPARENCY,' ! ',OB%REMOVABLE; CALL ADDSTR + WRITE(MYSTR,'(8I5,4F13.5,A,L1,1X,6I2)') OB%I1,OB%I2,OB%J1,OB%J2,OB%K1,OB%K2,OB%COLOR_INDICATOR,OB%TYPE_INDICATOR, & + REAL(OB%RGB,FB)/255._FB, OB%TRANSPARENCY,' ! ',OB%REMOVABLE,OB%EXPOSED_FACE_INDEX(1:6) ENDIF + CALL ADDSTR ENDDO ! Count circular vents diff --git a/Source/init.f90 b/Source/init.f90 index 981fc6861d..ea87bfeac2 100644 --- a/Source/init.f90 +++ b/Source/init.f90 @@ -771,7 +771,7 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM) DO J=OB%J1+1,OB%J2 I = OB%I1+1 ! Don't assign wall cell index to obstruction face pointing out of the computational domain - IF (I==1) CYCLE + IF (I==1) THEN ; OB%EXPOSED_FACE_INDEX(1)=1 ; CYCLE ; ENDIF IC = M%CELL_INDEX(I-1,J,K) IF (M%CELL(IC)%SOLID .AND. .NOT.M%OBSTRUCTION(M%CELL(IC)%OBST_INDEX)%REMOVABLE) CYCLE ! Permanently covered face IOR = -1 @@ -794,7 +794,7 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM) DO J=OB%J1+1,OB%J2 I = OB%I2 ! Don't assign wall cell index to obstruction face pointing out of the computational domain - IF (I==M%IBAR) CYCLE + IF (I==M%IBAR) THEN ; OB%EXPOSED_FACE_INDEX(2)=1 ; CYCLE ; ENDIF IC = M%CELL_INDEX(I+1,J,K) ! Permanently covered face IF (M%CELL(IC)%SOLID .AND. .NOT.M%OBSTRUCTION(M%CELL(IC)%OBST_INDEX)%REMOVABLE) CYCLE @@ -818,7 +818,7 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM) DO I=OB%I1+1,OB%I2 J = OB%J1+1 ! Don't assign wall cell index to obstruction face pointing out of the computational domain - IF (J==1) CYCLE + IF (J==1) THEN ; OB%EXPOSED_FACE_INDEX(3)=1 ; CYCLE ; ENDIF IC = M%CELL_INDEX(I,J-1,K) ! Permanently covered face IF (M%CELL(IC)%SOLID .AND. .NOT.M%OBSTRUCTION(M%CELL(IC)%OBST_INDEX)%REMOVABLE) CYCLE @@ -842,7 +842,7 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM) DO I=OB%I1+1,OB%I2 J = OB%J2 ! Don't assign wall cell index to obstruction face pointing out of the computational domain - IF (J==M%JBAR) CYCLE + IF (J==M%JBAR) THEN ; OB%EXPOSED_FACE_INDEX(4)=1 ; CYCLE ; ENDIF IC = M%CELL_INDEX(I,J+1,K) ! Permanently covered face IF (M%CELL(IC)%SOLID .AND. .NOT.M%OBSTRUCTION(M%CELL(IC)%OBST_INDEX)%REMOVABLE) CYCLE @@ -866,7 +866,7 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM) DO I=OB%I1+1,OB%I2 K = OB%K1+1 ! Don't assign wall cell index to obstruction face pointing out of the computational domain - IF (K==1) CYCLE + IF (K==1) THEN ; OB%EXPOSED_FACE_INDEX(5)=1 ; CYCLE ; ENDIF IC = M%CELL_INDEX(I,J,K-1) ! Permanently covered face IF (M%CELL(IC)%SOLID .AND. .NOT.M%OBSTRUCTION(M%CELL(IC)%OBST_INDEX)%REMOVABLE) CYCLE @@ -890,7 +890,7 @@ SUBROUTINE INITIALIZE_MESH_VARIABLES_1(DT,NM) DO I=OB%I1+1,OB%I2 K = OB%K2 ! Don't assign wall cell index to obstruction face pointing out of the computational domain - IF (K==M%KBAR) CYCLE + IF (K==M%KBAR) THEN ; OB%EXPOSED_FACE_INDEX(6)=1 ; CYCLE ; ENDIF IC = M%CELL_INDEX(I,J,K+1) ! Permanently covered face IF (M%CELL(IC)%SOLID .AND. .NOT.M%OBSTRUCTION(M%CELL(IC)%OBST_INDEX)%REMOVABLE) CYCLE @@ -2891,7 +2891,7 @@ SUBROUTINE INIT_WALL_CELL(NM,I,J,K,OBST_INDEX,IW,IOR,SURF_INDEX,IERR,TT) REAL(EB), INTENT(IN) :: TT REAL(EB) :: PX,PY,PZ,T_ACTIVATE,XIN,YIN,ZIN,DIST,XW,YW,ZW,RDN,AW,TSI,& ZZ_GET(1:N_TRACKED_SPECIES),RSUM_F,R1,RR,DELTA -INTEGER :: N,SURF_INDEX_NEW,IIG,JJG,KKG,IIO,JJO,KKO,IC,ICG,ICO,NOM_CHECK(0:1),BOUNDARY_TYPE +INTEGER :: N,SURF_INDEX_NEW,IIG,JJG,KKG,IIO,JJO,KKO,IC,ICG,ICO,NOM_CHECK(0:1),BOUNDARY_TYPE,FI LOGICAL :: VENT_FOUND,ALIGNED TYPE (MESH_TYPE), POINTER :: M,MM TYPE (OBSTRUCTION_TYPE), POINTER :: OBX @@ -3195,6 +3195,17 @@ SUBROUTINE INIT_WALL_CELL(NM,I,J,K,OBST_INDEX,IW,IOR,SURF_INDEX,IERR,TT) ENDIF CHECK_MESHES +! If this wall cell is attached to an OBST, check if the OBST face is exposed + +IF (OBST_INDEX>0) THEN + IF (.NOT.M%CELL(ICG)%SOLID .OR. M%OBSTRUCTION(M%CELL(ICG)%OBST_INDEX)%REMOVABLE) THEN + FI = ABS(IOR)*2 ; IF (IOR<0) FI = FI-1 + M%OBSTRUCTION(OBST_INDEX)%EXPOSED_FACE_INDEX(FI) = 1 + ENDIF +ENDIF + +! Ensure that the WALL_INDEX and SURF_INDEX can be identified from the abutting gas phase cell, ICG + M%CELL(ICG)%WALL_INDEX(-IOR) = IW M%CELL(ICG)%SURF_INDEX(-IOR) = SURF_INDEX_NEW diff --git a/Source/main.f90 b/Source/main.f90 index 16d0e60fbb..6c8ad37e70 100644 --- a/Source/main.f90 +++ b/Source/main.f90 @@ -78,23 +78,24 @@ PROGRAM FDS CALL OPENMP_INIT -! output version info if fds is invoked without any arguments -! (this must be done before MPI is initialized) +! Output version info if fds is invoked without any arguments. This must be done before MPI is initialized. CALL VERSION_INFO -! Initialize MPI (First executable lines of code) +! Initialize MPI CALL MPI_INIT_THREAD(REQUIRED,PROVIDED,IERR) CALL MPI_COMM_RANK(MPI_COMM_WORLD, MY_RANK, IERR) CALL MPI_COMM_SIZE(MPI_COMM_WORLD, N_MPI_PROCESSES, IERR) CALL MPI_GET_PROCESSOR_NAME(PNAME, PNAMELEN, IERR) +! Write out MPI process info to standard error (LU_ERR=0) + IF (MY_RANK==0) WRITE(LU_ERR,'(/A/)') ' Starting FDS ...' -CALL MPI_BARRIER(MPI_COMM_WORLD, IERR) +CALL MPI_BARRIER(MPI_COMM_WORLD,IERR) WRITE(LU_ERR,'(A,I6,A,A)') ' MPI Process ',MY_RANK,' started on ',PNAME(1:PNAMELEN) -CALL MPI_BARRIER(MPI_COMM_WORLD, IERR) +CALL MPI_BARRIER(MPI_COMM_WORLD,IERR) ! Check that MPI processes and OpenMP threads are working properly @@ -107,14 +108,13 @@ PROGRAM FDS CPU_TIME_START = CPUTIME ALLOCATE(T_USED(N_TIMERS)) ; T_USED = 0._EB ; T_USED(1) = CURRENT_TIME() -! Assign a compilation date (All Nodes) +! Assign a compilation date -CALL GET_INFO (REVISION,REVISION_DATE,COMPILE_DATE) +CALL GET_INFO(REVISION,REVISION_DATE,COMPILE_DATE) ! Read input from CHID.fds file and stop the code if any errors are found -CALL READ_DATA(DT) -CALL STOP_CHECK(1) +CALL READ_DATA(DT) ; CALL STOP_CHECK(1) IF (MY_RANK==0) THEN CALL WRITE_SUMMARY_INFO(LU_ERR,.TRUE.) @@ -137,27 +137,10 @@ PROGRAM FDS COUNTS,DISPLS,MPI_DOUBLE_PRECISION,MPI_COMM_WORLD,IERR) IF (MAXVAL(MAX_CELL_ASPECT_RATIO)>3.99_EB .AND. .NOT.CFL_VELOCITY_NORM_USER_SPECIFIED) CFL_VELOCITY_NORM=1 -! Open and write to Smokeview and status files +! Create output file names CALL ASSIGN_FILE_NAMES -CALL WRITE_SMOKEVIEW_FILE - -! Shut down the run if it is only for checking the set up - -IF (SETUP_ONLY .AND. .NOT.CHECK_MESH_ALIGNMENT) STOP_STATUS = SETUP_ONLY_STOP - -! Check for errors and shutdown if found - -CALL STOP_CHECK(1) - -! MPI process 0 reopens the Smokeview file for additional output - -IF (MY_RANK==0) THEN - OPEN(LU_SMV,FILE=FN_SMV,FORM='FORMATTED', STATUS='OLD',POSITION='APPEND') - CALL WRITE_STATUS_FILES -ENDIF - ! Start the clock T = T_BEGIN @@ -200,14 +183,25 @@ PROGRAM FDS IF (MY_RANK==0 .AND. VERBOSE) CALL VERBOSE_PRINTOUT('Completed INITIALIZE_MESH_VARIABLES_1') +! Write the Smokeview (.smv) file using parallel MPI writes + +CALL WRITE_SMOKEVIEW_FILE + ! Stop all the processes if this is just a set-up run -IF (CHECK_MESH_ALIGNMENT) THEN +IF (SETUP_ONLY .OR. CHECK_MESH_ALIGNMENT) THEN IF (MY_RANK==0) CALL INITIALIZE_DIAGNOSTIC_FILE(DT) STOP_STATUS = SETUP_ONLY_STOP IF (MY_RANK==0) WRITE(LU_ERR,'(A)') ' Checking mesh alignment. This could take a few tens of seconds...' + CALL STOP_CHECK(1) +ENDIF + +! MPI process 0 reopens the Smokeview file for additional output + +IF (MY_RANK==0) THEN + OPEN(LU_SMV,FILE=FN_SMV,FORM='FORMATTED', STATUS='OLD',POSITION='APPEND') + CALL WRITE_STATUS_FILES ENDIF -CALL STOP_CHECK(1) ! Allocate and initialize OMESH arrays to hold "other mesh" data for a given mesh diff --git a/Source/type.f90 b/Source/type.f90 index 5a471bc65a..91c47ace90 100644 --- a/Source/type.f90 +++ b/Source/type.f90 @@ -1051,9 +1051,10 @@ MODULE TYPES CHARACTER(LABEL_LENGTH) :: CTRL_ID='null' !< Name of controller CHARACTER(LABEL_LENGTH) :: ID='null' !< Name of obstruction - INTEGER, DIMENSION(-3:3) :: SURF_INDEX=0 !< SURFace properties for each face - INTEGER :: SURF_INDEX_INTERIOR=0 !< SURFace properties for a newly exposed interior cell - INTEGER, DIMENSION(3) :: RGB=(/0,0,0/) !< Color indices for Smokeview + INTEGER, DIMENSION(-3:3) :: SURF_INDEX=0 !< SURFace properties for each face + INTEGER, DIMENSION(1:6) :: EXPOSED_FACE_INDEX=0 !< OBST face exposed (1) or blocked (0) + INTEGER :: SURF_INDEX_INTERIOR=0 !< SURFace properties for a newly exposed interior cell + INTEGER, DIMENSION(3) :: RGB=(/0,0,0/) !< Color indices for Smokeview REAL(EB) :: TRANSPARENCY=1._EB !< Transparency index for Smokeview, 0=invisible, 1=solid REAL(EB) :: VOLUME_ADJUST=1._EB !< Effective volume divided by user specified volume