diff --git a/Manuals/FDS_User_Guide/FDS_User_Guide.tex b/Manuals/FDS_User_Guide/FDS_User_Guide.tex index a0777256589..7f2116abd61 100644 --- a/Manuals/FDS_User_Guide/FDS_User_Guide.tex +++ b/Manuals/FDS_User_Guide/FDS_User_Guide.tex @@ -10219,6 +10219,13 @@ \subsection{Heat Flux} The quantity {\ct RADIOMETER} is based on the {\em ellipsoidal radiometer} concept first proposed by Nils-Erik Gunners and described in the Ref.~\cite{Murthy:JResNIST2003}. The intent of such a device is to eliminate the contribution of convection from the measurement. +\item[\ct 'RADIOMETER GAS'] +\hspace{1in} \newline This output quantity is the same as {\ct 'RADIOMETER'} described above, except it can located anywhere within the computational domain and not just at a solid surface. It also has an arbitrary {\ct ORIENTATION} vector that points in any desired direction, much like a heat flux gauge. The {\ct ORIENTATION} vector need not be normalized, as in the following: +\begin{lstlisting} + &DEVC ID='hf', QUANTITY='RADIOMETER GAS', XYZ=..., ORIENTATION=-1,1,0/ +\end{lstlisting} +Note that the parameter {\ct SPATIAL\_STATISTIC} is not appropriate for this quantity, meaning that you cannot integrate this quantity over a plane or volume. However, you can use the parameter {\ct POINTS} to create a one-dimensional array of these devices (see Sec.~\ref{info:line_file}). + \item[\ct 'RADIATIVE HEAT FLUX GAS'] \hspace{1in} \newline This output quantity is the same as {\ct 'RADIATIVE HEAT FLUX'} described above, except it can located anywhere within the computational domain and not just at a solid surface. It also has an arbitrary {\ct ORIENTATION} vector that points in any desired direction, much like a heat flux gauge. The {\ct ORIENTATION} vector need not be normalized, as in the following: \begin{lstlisting} @@ -10239,6 +10246,8 @@ \subsection{Heat Flux} \end{labeling} +For {\ct QUANTITY} types ending in {\ct GAS}, a view angle for the radiation component of the device can be set using {\ct VIEW\_ANGLE} on {\PROP}. Wall devices can only have a 180$^\circ$ view angle. + \subsection{Adiabatic Surface Temperature} \label{info:AST} @@ -12622,6 +12631,7 @@ \section{\texorpdfstring{{\tt PROP}}{PROP} (Device Properties)} {\ct SPRAY\_PATTERN\_TABLE} & Character & Section~\ref{info:sprinklers} & & \\ \hline {\ct TIME\_CONSTANT} & Real & Section~\ref{info:THERMOCOUPLE} & s & \\ \hline {\ct VELOCITY\_COMPONENT} & Integer & Section~\ref{info:velocity_patch} & & \\ \hline +{\ct VIEW\_ANGLE} & Real & Section~\ref{info:heat_flux} & degrees & 180. \\ \hline \end{longtable} diff --git a/Source/cons.f90 b/Source/cons.f90 index e6199b732bf..63397344a75 100644 --- a/Source/cons.f90 +++ b/Source/cons.f90 @@ -714,7 +714,9 @@ MODULE GLOBAL_CONSTANTS REAL(EB), POINTER, DIMENSION(:,:) :: ORIENTATION_VECTOR !< Global array of orientation vectors INTEGER, ALLOCATABLE, DIMENSION(:) :: NEAREST_RADIATION_ANGLE !< Index of the rad angle most opposite the given ORIENTATION_VECTOR -INTEGER :: N_ORIENTATION_VECTOR !< Number of ORIENTATION_VECTORs +REAL(EB), POINTER, DIMENSION(:) :: ORIENTATION_VIEW_ANGLE !< View angle of the given ORIENTATION_VECTOR +REAL(EB), ALLOCATABLE, DIMENSION(:) :: VIEW_ANGLE_AREA !< View angle area ORIENTATION_VECTOR +INTEGER :: N_ORIENTATION_VECTOR !< Number of ORIENTATION_VECTORs INTEGER :: TGA_SURF_INDEX=-100 !< Surface properties to use for special TGA calculation INTEGER :: TGA_WALL_INDEX=-100 !< Wall index to use for special TGA calculation @@ -791,6 +793,7 @@ MODULE RADCONS IMPLICIT NONE (TYPE,EXTERNAL) REAL(EB), ALLOCATABLE, DIMENSION(:,:) :: DLN !< Wall-normal matrix +REAL(EB), ALLOCATABLE, DIMENSION(:,:) :: DLANG !< Angles REAL(EB), ALLOCATABLE, DIMENSION(:,:) :: ORIENTATION_FACTOR !< Fraction of radiation angle corresponding to a particular direction REAL(EB), ALLOCATABLE, DIMENSION(:) :: BBFRAC !< Fraction of blackbody radiation REAL(EB), ALLOCATABLE, DIMENSION(:) :: WL_LOW !< Lower wavelength limit of the spectral band diff --git a/Source/data.f90 b/Source/data.f90 index c6f82635aa2..bd518da5d37 100644 --- a/Source/data.f90 +++ b/Source/data.f90 @@ -1392,6 +1392,7 @@ SUBROUTINE DEFINE_OUTPUT_QUANTITIES OUTPUT_QUANTITY(-21)%SHORT_NAME= 'h' OUTPUT_QUANTITY(-22)%NAME= 'RADIOMETER' +OUTPUT_QUANTITY(-22)%OLD_NAME= 'RADIOMETER GAS' OUTPUT_QUANTITY(-22)%UNITS= 'kW/m2' OUTPUT_QUANTITY(-22)%SHORT_NAME= 'radio' diff --git a/Source/devc.f90 b/Source/devc.f90 index 4deb36f9170..35b9d9fad6d 100644 --- a/Source/devc.f90 +++ b/Source/devc.f90 @@ -14,7 +14,7 @@ MODULE DEVICE_VARIABLES ACTIVATION_TEMPERATURE,ACTIVATION_OBSCURATION, & ALPHA_E,ALPHA_C,BETA_E,BETA_C,CHARACTERISTIC_VELOCITY,PARTICLE_VELOCITY,MASS_FLOW_RATE,FLOW_RATE,FLOW_TAU, & GAUGE_EMISSIVITY,GAUGE_TEMPERATURE,INITIAL_TEMPERATURE,K_FACTOR,C_FACTOR,OPERATING_PRESSURE,OFFSET,& - SPRAY_ANGLE(2,2),P0=0._EB,PX(3)=0._EB,PXX(3,3)=0._EB + SPRAY_ANGLE(2,2),P0=0._EB,PX(3)=0._EB,PXX(3,3)=0._EB,VIEW_ANGLE INTEGER :: PDPA_M=0,PDPA_N=0,N_SMOKEVIEW_PARAMETERS=0,N_SMOKEVIEW_IDS=0,N_INSERT,I_VEL=0,PARTICLES_PER_SECOND LOGICAL :: PDPA_INTEGRATE=.TRUE.,PDPA_NORMALIZE=.TRUE.,HISTOGRAM_NORMALIZE=.TRUE.,HISTOGRAM=.FALSE., & HISTOGRAM_CUMULATIVE=.FALSE.,SPARK=.FALSE. diff --git a/Source/dump.f90 b/Source/dump.f90 index cc513f72739..9aa7069c912 100644 --- a/Source/dump.f90 +++ b/Source/dump.f90 @@ -3660,11 +3660,11 @@ SUBROUTINE WRITE_DIAGNOSTICS(T,DT) CALL GET_DATE_ISO_8601(DATE) CALL CPU_TIME(CPUTIME) IF (ABS(T)<=999.99999_EB) THEN - WRITE(LU_STEPS,'(I7,",",A,",",E10.3,",",F10.5,",",E12.5)') ICYC,TRIM(DATE),DT,T,CPUTIME - CPU_TIME_START + WRITE(LU_STEPS,'(I8,",",A,",",E10.3,",",F10.5,",",E12.5)') ICYC,TRIM(DATE),DT,T,CPUTIME - CPU_TIME_START ELSEIF (ABS(T)>999.99999_EB .AND. ABS(T)<=99999.999_EB) THEN - WRITE(LU_STEPS,'(I7,",",A,",",E10.3,",",F10.3,",",E12.5)') ICYC,TRIM(DATE),DT,T,CPUTIME - CPU_TIME_START + WRITE(LU_STEPS,'(I8,",",A,",",E10.3,",",F10.3,",",E12.5)') ICYC,TRIM(DATE),DT,T,CPUTIME - CPU_TIME_START ELSE - WRITE(LU_STEPS,'(I7,",",A,",",E10.3,",",F10.1,",",E12.5)') ICYC,TRIM(DATE),DT,T,CPUTIME - CPU_TIME_START + WRITE(LU_STEPS,'(I8,",",A,",",E10.3,",",F10.1,",",E12.5)') ICYC,TRIM(DATE),DT,T,CPUTIME - CPU_TIME_START ENDIF ! Write abridged output to the .err file @@ -3672,26 +3672,26 @@ SUBROUTINE WRITE_DIAGNOSTICS(T,DT) IF (ABS(TIME_SHRINK_FACTOR-1._EB) < TWO_EPSILON_EB) THEN IF (ABS(T)<=0.0001) THEN - WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I7,A,F10.5,A)') 'Time Step:',ICYC,', Simulation Time:',T,' s' + WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I8,A,F10.5,A)') 'Time Step:',ICYC,', Simulation Time:',T,' s' ELSEIF (ABS(T)>0.0001 .AND. ABS(T) <=0.001) THEN - WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I7,A,F10.4,A)') 'Time Step:',ICYC,', Simulation Time:',T,' s' + WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I8,A,F10.4,A)') 'Time Step:',ICYC,', Simulation Time:',T,' s' ELSEIF (ABS(T)>0.001 .AND. ABS(T)<=0.01) THEN - WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I7,A,F10.3,A)') 'Time Step:',ICYC,', Simulation Time:',T,' s' + WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I8,A,F10.3,A)') 'Time Step:',ICYC,', Simulation Time:',T,' s' ELSE - WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I7,A,F10.2,A)') 'Time Step:',ICYC,', Simulation Time:',T,' s' + WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I8,A,F10.2,A)') 'Time Step:',ICYC,', Simulation Time:',T,' s' ENDIF ELSE STIME = T_BEGIN + (T-T_BEGIN) * TIME_SHRINK_FACTOR DTS = DT * TIME_SHRINK_FACTOR IF (ABS(STIME)<=0.0001) THEN - WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I7,A,F10.5,A)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,' s' + WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I8,A,F10.5,A)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,' s' ELSEIF (ABS(STIME)>0.0001 .AND. ABS(STIME) <=0.001) THEN - WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I7,A,F10.4,A)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,' s' + WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I8,A,F10.4,A)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,' s' ELSEIF (ABS(STIME)>0.001 .AND. ABS(STIME)<=0.01) THEN - WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I7,A,F10.3,A)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,' s' + WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I8,A,F10.3,A)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,' s' ELSE - WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I7,A,F10.2,A)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,' s' + WRITE(SIMPLE_OUTPUT_ERR,'(1X,A,I8,A,F10.2,A)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,' s' ENDIF ENDIF @@ -3711,19 +3711,19 @@ SUBROUTINE WRITE_DIAGNOSTICS(T,DT) IF (ABS(TIME_SHRINK_FACTOR-1._EB) < TWO_EPSILON_EB) THEN IF (ABS(T)<=0.0001) THEN - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.6,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.6,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& ' s, Pressure Iterations: ',PRESSURE_ITERATIONS ELSEIF (ABS(T)>0.0001 .AND. ABS(T) <=0.001) THEN - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.5,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.5,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& ' s, Pressure Iterations: ',PRESSURE_ITERATIONS ELSEIF (ABS(T)>0.001 .AND. ABS(T)<=0.01) THEN - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.4,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.4,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& ' s, Pressure Iterations: ',PRESSURE_ITERATIONS ELSEIF (ABS(T)>0.01 .AND. ABS(T)<=0.1) THEN - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.3,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.3,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& ' s, Pressure Iterations: ',PRESSURE_ITERATIONS ELSE - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.2,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.2,A,F8.5,A,I0)') 'Time Step:',ICYC,', Simulation Time:',T,' s, Step Size:',DT,& ' s, Pressure Iterations: ',PRESSURE_ITERATIONS ENDIF @@ -3732,19 +3732,19 @@ SUBROUTINE WRITE_DIAGNOSTICS(T,DT) ELSE IF (ABS(STIME)<=0.0001) THEN - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.6,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.6,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& ' s, Scaled Step Size:',DTS,' s, Pressure Iterations: ',PRESSURE_ITERATIONS ELSEIF (ABS(STIME)>0.0001 .AND. ABS(STIME) <=0.001) THEN - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.5,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.5,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& ' s, Scaled Step Size:',DTS,' s, Pressure Iterations: ',PRESSURE_ITERATIONS ELSEIF (ABS(STIME)>0.001 .AND. ABS(STIME)<=0.01) THEN - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.4,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.4,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& ' s, Scaled Step Size:',DTS,' s, Pressure Iterations: ',PRESSURE_ITERATIONS ELSEIF (ABS(STIME)>0.01 .AND. ABS(STIME)<=0.1) THEN - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.3,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.3,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& ' s, Scaled Step Size:',DTS,' s, Pressure Iterations: ',PRESSURE_ITERATIONS ELSE - WRITE(SIMPLE_OUTPUT,'(1X,A,I7,A,F10.2,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& + WRITE(SIMPLE_OUTPUT,'(1X,A,I8,A,F10.2,A,F8.5,A,I0)') 'Time Step:',ICYC,', Scaled Simulation Time:',STIME,& ' s, Scaled Step Size:',DTS,' s, Pressure Iterations: ',PRESSURE_ITERATIONS ENDIF ENDIF @@ -3754,7 +3754,7 @@ SUBROUTINE WRITE_DIAGNOSTICS(T,DT) ! Detailed diagnostics to the .out file CALL GET_DATE(DATE) -WRITE(LU_OUTPUT,'(7X,A,I7,3X,A)') 'Time Step ',ICYC,TRIM(DATE) +WRITE(LU_OUTPUT,'(7X,A,I8,3X,A)') 'Time Step ',ICYC,TRIM(DATE) IF (ABS(TIME_SHRINK_FACTOR-1._EB) < TWO_EPSILON_EB) THEN IF (ABS(T)<=0.0001) THEN WRITE(LU_OUTPUT,150) DT,T diff --git a/Source/part.f90 b/Source/part.f90 index 5764464b5f6..3bddb338505 100644 --- a/Source/part.f90 +++ b/Source/part.f90 @@ -1500,6 +1500,7 @@ SUBROUTINE VOLUME_INIT_PARTICLE IF (DV%QUANTITY(1)=='RADIATIVE HEAT FLUX' .OR. & DV%QUANTITY(1)=='GAUGE HEAT FLUX' .OR. & DV%QUANTITY(1)=='RADIANCE' .OR. & + DV%QUANTITY(1)=='RADIOMETER' .OR. & DV%QUANTITY(1)=='ADIABATIC SURFACE TEMPERATURE') THEN IF (LPC%ID=='RESERVED TARGET PARTICLE') THEN ! use the orientation of the DEVC LP%ORIENTATION_INDEX = DV%ORIENTATION_INDEX diff --git a/Source/radi.f90 b/Source/radi.f90 index 70e79868023..6b05a74987b 100644 --- a/Source/radi.f90 +++ b/Source/radi.f90 @@ -2787,7 +2787,8 @@ SUBROUTINE INIT_RADIATION USE MIEV USE RADCAL_CALC USE WSGG_ARRAYS -REAL(EB) :: THETAUP,THETALOW,PHIUP,PHILOW,F_THETA,PLANCK_C2,KSI,LT,RCRHO,YY,YY2,BBF,AP0,AMEAN,RADIANCE,TRANSMISSIVITY,X_N2 +REAL(EB) :: THETAUP,THETALOW,PHIUP,PHILOW,F_THETA,PLANCK_C2,KSI,LT,RCRHO,YY,YY2,BBF,AP0,AMEAN,RADIANCE,TRANSMISSIVITY,X_N2,& + THETA,PHI INTEGER :: N,I,J,K,IPC,IZERO,NN,NI,II,JJ,IIM,JJM,IBND,NS,NS2,NRA,NSB,RADCAL_TEMP(16)=0,RCT_SKIP=-1,OR_IN,I1,I2,IO TYPE (LAGRANGIAN_PARTICLE_CLASS_TYPE), POINTER :: LPC REAL(EB), ALLOCATABLE, DIMENSION(:) :: COSINE_ARRAY @@ -2821,6 +2822,8 @@ SUBROUTINE INIT_RADIATION CALL ChkMemErr('RADI','DLN',IZERO) ALLOCATE(DLM(1:NRA,3),STAT=IZERO) CALL ChkMemErr('RADI','DLM',IZERO) +ALLOCATE(DLANG(3,1:NRA),STAT=IZERO) +CALL ChkMemErr('RADI','DLANG',IZERO) ! Determine mean direction normals and sweeping orders ! as described in the FDS Tech. Ref. Guide Vol. 1 Sec. 6.2.2. @@ -2832,6 +2835,7 @@ SUBROUTINE INIT_RADIATION THETALOW = PI*REAL(I-1)/REAL(NRT) THETAUP = PI*REAL(I)/REAL(NRT) F_THETA = 0.5_EB*(THETAUP-THETALOW - COS(THETAUP)*SIN(THETAUP) + COS(THETALOW)*SIN(THETALOW)) + THETA = 0.5_EB*(THETAUP+THETALOW) IF (CYLINDRICAL) THEN PHILOW = PI*REAL(J-1)/REAL(NRP(I)) PHIUP = PI*REAL(J)/REAL(NRP(I)) @@ -2842,21 +2846,31 @@ SUBROUTINE INIT_RADIATION PHILOW = TWOPI*REAL(J-1)/REAL(NRP(I)) PHIUP = TWOPI*REAL(J)/REAL(NRP(I)) ENDIF + PHI=0.5_EB*(PHILOW+PHIUP) RSA(N) = (PHIUP-PHILOW)*(COS(THETALOW)-COS(THETAUP)) IF (CYLINDRICAL) THEN DLX(N) = (SIN(PHIUP)-SIN(PHILOW)) *F_THETA DLY(N) = (-SIN(DPHI0/2.)*(SIN(PHIUP)-SIN(PHILOW)) +COS(DPHI0/2.)*(COS(PHILOW)-COS(PHIUP)))*F_THETA DLB(N) = (-SIN(DPHI0/2.)*(SIN(PHIUP)-SIN(PHILOW)) -COS(DPHI0/2.)*(COS(PHILOW)-COS(PHIUP)))*F_THETA DLZ(N) = 0.5_EB*(PHIUP-PHILOW) * ((SIN(THETAUP))**2-(SIN(THETALOW))**2) + DLANG(1,N) = SIN(THETA)*COS(PHI) + DLANG(2,N) = SIN(THETA)*SIN(PHI) + DLANG(3,N) = COS(THETA) IF (N==1000000) WRITE(LU_ERR,'(A)') 'This line should never get executed. It is here only to prevent optimization.' ELSEIF (TWO_D) THEN DLX(N) = (SIN(PHIUP)-SIN(PHILOW))*F_THETA DLY(N) = 0._EB DLZ(N) = (COS(PHILOW)-COS(PHIUP))*F_THETA + DLANG(1,N) = COS(PHI) + DLANG(2,N) = 0._EB + DLANG(3,N) = SIN(PHI) ELSE DLX(N) = (SIN(PHIUP)-SIN(PHILOW))*F_THETA DLY(N) = (COS(PHILOW)-COS(PHIUP))*F_THETA DLZ(N) = 0.5_EB*(PHIUP-PHILOW) * ((SIN(THETAUP))**2-(SIN(THETALOW))**2) + DLANG(1,N) = SIN(THETA)*COS(PHI) + DLANG(2,N) = SIN(THETA)*SIN(PHI) + DLANG(3,N) = COS(THETA) ENDIF ENDDO ENDDO @@ -3392,14 +3406,22 @@ SUBROUTINE INIT_RADIATION ALLOCATE(COSINE_ARRAY(1:NRA)) ALLOCATE(NEAREST_RADIATION_ANGLE(N_ORIENTATION_VECTOR)) + ALLOCATE(VIEW_ANGLE_AREA(N_ORIENTATION_VECTOR)) + VIEW_ANGLE_AREA = 0._EB DO IO=1,N_ORIENTATION_VECTOR DO N=1,NRA COSINE_ARRAY(N) = ORIENTATION_VECTOR(1,IO)*DLX(N) + & ORIENTATION_VECTOR(2,IO)*DLY(N) + & ORIENTATION_VECTOR(3,IO)*DLZ(N) + IF (-(ORIENTATION_VECTOR(1,IO)*DLANG(1,N) + & + ORIENTATION_VECTOR(2,IO)*DLANG(2,N) + & + ORIENTATION_VECTOR(3,IO)*DLANG(3,N)) > ORIENTATION_VIEW_ANGLE(IO)) & + VIEW_ANGLE_AREA(IO) = VIEW_ANGLE_AREA(IO) - COSINE_ARRAY(N) ENDDO NEAREST_RADIATION_ANGLE(IO) = MINLOC(COSINE_ARRAY,DIM=1) + VIEW_ANGLE_AREA(IO) = PI/VIEW_ANGLE_AREA(IO) ENDDO + DEALLOCATE(COSINE_ARRAY) ENDIF @@ -4391,18 +4413,19 @@ SUBROUTINE RADIATION_FVM TEMP_ORIENTATION = TEMP_ORIENTATION / & (SQRT(TEMP_ORIENTATION(1)**2+TEMP_ORIENTATION(2)**2+TEMP_ORIENTATION(3)**2) & +TWO_EPSILON_EB) - COS_DL = -(TEMP_ORIENTATION(1)*DLX(N) + & - TEMP_ORIENTATION(2)*DLY(N) + & - TEMP_ORIENTATION(3)*DLZ(N)) - IF (COS_DL>0._EB) THEN + COS_DL = -DOT_PRODUCT(TEMP_ORIENTATION(1:3),DLANG(1:3,N)) + IF (COS_DL>ORIENTATION_VIEW_ANGLE(LP%ORIENTATION_INDEX)) THEN + COS_DL = -(TEMP_ORIENTATION(1)*DLX(N) + & + TEMP_ORIENTATION(2)*DLY(N) + & + TEMP_ORIENTATION(3)*DLZ(N)) BR => BOUNDARY_RADIA(LP%BR_INDEX) IF (LPC%MASSLESS_TARGET) THEN - BR%BAND(IBND)%ILW(N) = COS_DL * IL(BC%IIG,BC%JJG,BC%KKG) + BR%BAND(IBND)%ILW(N) = COS_DL * IL(BC%IIG,BC%JJG,BC%KKG) * VIEW_ANGLE_AREA(LP%ORIENTATION_INDEX) IF (N==NEAREST_RADIATION_ANGLE(LP%ORIENTATION_INDEX)) & BR%IL(IBND) = IL(BC%IIG,BC%JJG,BC%KKG) ELSE ! IL_UP does not account for the absorption of radiation within the cell occupied by the particle - BR%BAND(IBND)%ILW(N) = COS_DL * IL_UP(BC%IIG,BC%JJG,BC%KKG) + BR%BAND(IBND)%ILW(N) = COS_DL * IL_UP(BC%IIG,BC%JJG,BC%KKG) * VIEW_ANGLE_AREA(LP%ORIENTATION_INDEX) ENDIF ENDIF CYCLE PARTICLE_RADIATION_LOOP @@ -4412,18 +4435,19 @@ SUBROUTINE RADIATION_FVM CASE(0) CYCLE PARTICLE_RADIATION_LOOP CASE(1) - COS_DL = -(ORIENTATION_VECTOR(1,LP%ORIENTATION_INDEX)*DLX(N) + & - ORIENTATION_VECTOR(2,LP%ORIENTATION_INDEX)*DLY(N) + & - ORIENTATION_VECTOR(3,LP%ORIENTATION_INDEX)*DLZ(N)) - IF (COS_DL>0._EB) THEN + COS_DL = -DOT_PRODUCT(ORIENTATION_VECTOR(1:3,LP%ORIENTATION_INDEX),DLANG(1:3,N)) + IF (COS_DL>ORIENTATION_VIEW_ANGLE(LP%ORIENTATION_INDEX)) THEN + COS_DL = -(ORIENTATION_VECTOR(1,LP%ORIENTATION_INDEX)*DLX(N) + & + ORIENTATION_VECTOR(2,LP%ORIENTATION_INDEX)*DLY(N) + & + ORIENTATION_VECTOR(3,LP%ORIENTATION_INDEX)*DLZ(N)) BR => BOUNDARY_RADIA(LP%BR_INDEX) IF (LPC%MASSLESS_TARGET) THEN - BR%BAND(IBND)%ILW(N) = COS_DL * IL(BC%IIG,BC%JJG,BC%KKG) + BR%BAND(IBND)%ILW(N) = COS_DL * IL(BC%IIG,BC%JJG,BC%KKG) * VIEW_ANGLE_AREA(LP%ORIENTATION_INDEX) IF (N==NEAREST_RADIATION_ANGLE(LP%ORIENTATION_INDEX)) & BR%IL(IBND) = IL(BC%IIG,BC%JJG,BC%KKG) ELSE ! IL_UP does not account for the absorption of radiation within the cell occupied by the particle - BR%BAND(IBND)%ILW(N) = COS_DL * IL_UP(BC%IIG,BC%JJG,BC%KKG) + BR%BAND(IBND)%ILW(N) = COS_DL * IL_UP(BC%IIG,BC%JJG,BC%KKG) * VIEW_ANGLE_AREA(LP%ORIENTATION_INDEX) ENDIF ENDIF END SELECT diff --git a/Source/read.f90 b/Source/read.f90 index 532cda91955..bd69c8dae81 100644 --- a/Source/read.f90 +++ b/Source/read.f90 @@ -13,7 +13,7 @@ MODULE READ_INPUT USE MESH_POINTERS USE OUTPUT_DATA USE COMP_FUNCTIONS, ONLY: CHECKREAD, SHUTDOWN, CHECK_XB, SCAN_INPUT_FILE -USE MEMORY_FUNCTIONS, ONLY: ChkMemErr,REALLOCATE2D +USE MEMORY_FUNCTIONS, ONLY: ChkMemErr,REALLOCATE2D,REALLOCATE USE COMP_FUNCTIONS, ONLY: GET_INPUT_FILE USE MISC_FUNCTIONS, ONLY: SEARCH_CONTROLLER,WRITE_SUMMARY_INFO USE HVAC_ROUTINES, ONLY: READ_HVAC,PROC_HVAC @@ -96,7 +96,9 @@ SUBROUTINE READ_DATA(DT) N_ORIENTATION_VECTOR = 0 ALLOCATE(ORIENTATION_VECTOR(3,0:10)) +ALLOCATE(ORIENTATION_VIEW_ANGLE(0:10)) ORIENTATION_VECTOR(1:3,0) = (/0._EB,0._EB,-1._EB/) +ORIENTATION_VIEW_ANGLE = 0._EB ! Open the input file @@ -5319,7 +5321,7 @@ SUBROUTINE READ_PART N_LAGRANGIAN_CLASSES = N_LAGRANGIAN_CLASSES_READ -! Add reserved INIT lines to account for devices for 'RADIATIVE HEAT FLUX GAS' or 'ADIABATIC SURFACE TEMPERATURE GAS' +! Add reserved INIT lines to account for devices for 'RADIATIVE HEAT FLUX GAS', 'ADIABATIC SURFACE TEMPERATURE GAS', etc. IF (TARGET_PARTICLES_INCLUDED) N_LAGRANGIAN_CLASSES = N_LAGRANGIAN_CLASSES + 1 @@ -5633,8 +5635,10 @@ SUBROUTINE READ_PART LPC%ORIENTATION_INDEX = N_ORIENTATION_VECTOR IF (N_ORIENTATION_VECTOR>UBOUND(ORIENTATION_VECTOR,DIM=2)) THEN ORIENTATION_VECTOR => REALLOCATE2D(ORIENTATION_VECTOR,1,3,0,N_ORIENTATION_VECTOR+10) + ORIENTATION_VIEW_ANGLE => REALLOCATE(ORIENTATION_VIEW_ANGLE,0,N_ORIENTATION_VECTOR+10) ENDIF ORIENTATION_VECTOR(1:3,N_ORIENTATION_VECTOR) = ORIENTATION(1:3)/ NORM2(ORIENTATION) + ORIENTATION_VIEW_ANGLE(N_ORIENTATION_VECTOR) = 0._EB ENDIF LPC%FREE_AREA_FRACTION = FREE_AREA_FRACTION LPC%POROUS_VOLUME_FRACTION = POROUS_VOLUME_FRACTION @@ -5973,7 +5977,7 @@ SUBROUTINE READ_PROP C_FACTOR,CHARACTERISTIC_VELOCITY,ORIFICE_DIAMETER,EMISSIVITY, & PARTICLE_VELOCITY,FLOW_RATE,FLOW_TAU,GAUGE_EMISSIVITY,GAUGE_TEMPERATURE,INITIAL_TEMPERATURE,K_FACTOR,& LENGTH,SPRAY_ANGLE(2,2),OFFSET,OPERATING_PRESSURE,RTI,PDPA_START,PDPA_END,PDPA_RADIUS,MASS_FLOW_RATE,& - SPRAY_PATTERN_MU,SPRAY_PATTERN_BETA,HISTOGRAM_LIMITS(2),P0,PX(3),PXX(3,3),TIME_CONSTANT + SPRAY_PATTERN_MU,SPRAY_PATTERN_BETA,HISTOGRAM_LIMITS(2),P0,PX(3),PXX(3,3),TIME_CONSTANT,VIEW_ANGLE INTEGER ::I,N,NN,PDPA_M,PDPA_N,PARTICLES_PER_SECOND,VELOCITY_COMPONENT,HISTOGRAM_NBINS,FED_ACTIVITY LOGICAL :: PDPA_INTEGRATE,PDPA_NORMALIZE,HISTOGRAM_NORMALIZE,HISTOGRAM,HISTOGRAM_CUMULATIVE,SPARK CHARACTER(LABEL_LENGTH) :: SMOKEVIEW_ID(SMOKEVIEW_OBJECTS_DIMENSION),QUANTITY='null',PART_ID='null',FLOW_RAMP='null', & @@ -5991,7 +5995,8 @@ SUBROUTINE READ_PROP PDPA_INTEGRATE,PDPA_M,PDPA_N,PDPA_NORMALIZE,PDPA_RADIUS,& PDPA_START,PRESSURE_RAMP,PX,PXX,QUANTITY,RTI,SMOKEVIEW_ID,SMOKEVIEW_PARAMETERS,SPARK,& SPEC_ID,SPECIFIC_HEAT,SPRAY_ANGLE,& - SPRAY_PATTERN_BETA,SPRAY_PATTERN_MU,SPRAY_PATTERN_SHAPE,SPRAY_PATTERN_TABLE,TIME_CONSTANT,VELOCITY_COMPONENT + SPRAY_PATTERN_BETA,SPRAY_PATTERN_MU,SPRAY_PATTERN_SHAPE,SPRAY_PATTERN_TABLE,TIME_CONSTANT,VELOCITY_COMPONENT,& + VIEW_ANGLE ! Count the PROP lines in the input file. Note how many of these are cables. @@ -6255,6 +6260,11 @@ SUBROUTINE READ_PROP CALL SHUTDOWN(MESSAGE) ; RETURN ENDIF ENDIF + PY%VIEW_ANGLE = VIEW_ANGLE + IF (PY%VIEW_ANGLE180._EB) THEN + WRITE(MESSAGE,'(A,A,A)') 'ERROR(xxx): PROP ',TRIM(PY%ID),' VIEW_ANGLE must be between 0 and 180.' + CALL SHUTDOWN(MESSAGE) ; RETURN + ENDIF ENDDO READ_PROP_LOOP @@ -6328,6 +6338,7 @@ SUBROUTINE SET_PROP_DEFAULTS TIME_CONSTANT = -1._EB FED_ACTIVITY = 2 ! light work VELOCITY_COMPONENT = 0 +VIEW_ANGLE = 180._EB END SUBROUTINE SET_PROP_DEFAULTS END SUBROUTINE READ_PROP @@ -11899,7 +11910,7 @@ SUBROUTINE READ_INIT ELSE - ! Use information from DEVC line to create an INIT line for 'RADIATIVE HEAT FLUX GAS' or 'ADIABATIC SURFACE TEMPERATURE GAS' + ! Use information from DEVC line to create an INIT line for 'RADIATIVE HEAT FLUX GAS', 'ADIABATIC SURFACE TEMPERATURE GAS', etc. CALL SET_INIT_DEFAULTS DV => DEVICE(INIT_RESERVED(N-N_INIT_READ)%DEVC_INDEX) ! First device in the line of POINTS @@ -12595,7 +12606,7 @@ SUBROUTINE READ_DEVC REAL(EB) :: DEPTH,ORIENTATION(3),ROTATION,SETPOINT,FLOWRATE,BYPASS_FLOWRATE,DELAY,XYZ(3),CONVERSION_FACTOR,CONVERSION_ADDEND, & SMOOTHING_FACTOR,QUANTITY_RANGE(2),STATISTICS_START,STATISTICS_END,COORD_FACTOR,CELL_L,& TIME_PERIOD,FORCE_DIRECTION(3),XI,YJ,ZK,XBP(6),DX,DY,DZ,& - POINTS_ARRAY_X(POINTS_ARRAY_DIM),POINTS_ARRAY_Y(POINTS_ARRAY_DIM),POINTS_ARRAY_Z(POINTS_ARRAY_DIM) + POINTS_ARRAY_X(POINTS_ARRAY_DIM),POINTS_ARRAY_Y(POINTS_ARRAY_DIM),POINTS_ARRAY_Z(POINTS_ARRAY_DIM),VIEW_ANGLE CHARACTER(LABEL_LENGTH) :: QUANTITY,QUANTITY2,PROP_ID,CTRL_ID,DEVC_ID,INIT_ID,SURF_ID,SPATIAL_STATISTIC,TEMPORAL_STATISTIC,& MOVE_ID,STATISTICS,PART_ID,MATL_ID,SPEC_ID,UNITS, & DUCT_ID,NODE_ID(2),D_ID,R_ID,X_ID,Y_ID,Z_ID,NO_UPDATE_DEVC_ID,NO_UPDATE_CTRL_ID,REAC_ID,XYZ_UNITS @@ -12644,7 +12655,7 @@ SUBROUTINE READ_DEVC ALLOCATE(DEVICE(N_DEVC),STAT=IZERO) ; CALL ChkMemErr('READ','DEVICE',IZERO) -! Speceial case for QUANTITY='RADIATIVE HEAT FLUX GAS' or 'ADIABATIC SURFACE TEMPERATURE GAS' +! Speceial case for QUANTITY='RADIATIVE HEAT FLUX GAS', 'ADIABATIC SURFACE TEMPERATURE GAS', etc. ALLOCATE(INIT_RESERVED(N_DEVC),STAT=IZERO) ; CALL ChkMemErr('READ','INIT_RESERVED',IZERO) N_INIT_RESERVED = 0 @@ -12809,7 +12820,8 @@ SUBROUTINE READ_DEVC DO I=1,N_ORIENTATION_VECTOR IF (ABS(ORIENTATION(1)-ORIENTATION_VECTOR(1,I))UBOUND(ORIENTATION_VECTOR,DIM=2)) THEN ORIENTATION_VECTOR => REALLOCATE2D(ORIENTATION_VECTOR,1,3,0,N_ORIENTATION_VECTOR+10) + ORIENTATION_VIEW_ANGLE => REALLOCATE(ORIENTATION_VIEW_ANGLE,0,N_ORIENTATION_VECTOR+10) ENDIF IF (ALL(ABS(ORIENTATION(1:3)) 0) THEN + IF (PROPERTY(DV%PROP_INDEX)%VIEW_ANGLE < 180._EB) & + ORIENTATION_VIEW_ANGLE(DV%ORIENTATION_INDEX) = COS(PROPERTY(DV%PROP_INDEX)%VIEW_ANGLE/360._EB * PI) + ENDIF ! Create an auto-ignition exclusion zone (AIT) in the cell containing a SPARK IF (PROPERTY(DV%PROP_INDEX)%SPARK .AND. PROCESS(DV%MESH)==MY_RANK) THEN