Merge pull request #13481 from mcgratta/master

FDS Source: Issue #13440. Clarify some radiation arrays

Source/cons.f90

@@ -727,8 +727,8 @@ 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
-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
+REAL(EB), POINTER, DIMENSION(:) :: COS_HALF_VIEW_ANGLE     !< View angle of the given ORIENTATION_VECTOR
+REAL(EB), ALLOCATABLE, DIMENSION(:) :: VIEW_ANGLE_FACTOR        !< View angle area ORIENTATION_VECTOR
 INTEGER :: N_ORIENTATION_VECTOR                               !< Number of ORIENTATION_VECTORs
 
 INTEGER :: TGA_MESH_INDEX=HUGE(INTEGER_ONE)  !< Mesh for the special TGA calculation

Source/radi.f90

@@ -2788,7 +2788,7 @@ SUBROUTINE INIT_RADIATION
 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,&
-            THETA,PHI
+            THETA,PHI,DLO
 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,IO
 TYPE (LAGRANGIAN_PARTICLE_CLASS_TYPE), POINTER :: LPC
 REAL(EB), ALLOCATABLE, DIMENSION(:) :: COSINE_ARRAY
@@ -3389,24 +3389,26 @@ SUBROUTINE INIT_RADIATION
 ENDDO
 
 ! Determine angle factors for Lagrangian particles with ORIENTATION
+! COSINE_ARRAY holds the cosines of the angles formed by the orientation vector and the radiation directions.
+! DLO is the integral of the orientation vector dotted with the directional solid angle of the radiation directions.
+! VIEW_ANGLE_FACTOR is the reduction of the radiation due to a view angle less than 180, like a narrow field of view radiometer.
 
 IF (SOLID_PARTICLES) THEN
    ALLOCATE(COSINE_ARRAY(1:NRA))
    ALLOCATE(NEAREST_RADIATION_ANGLE(N_ORIENTATION_VECTOR))
-   ALLOCATE(VIEW_ANGLE_AREA(N_ORIENTATION_VECTOR))
-   VIEW_ANGLE_AREA = 0._EB
+   ALLOCATE(VIEW_ANGLE_FACTOR(N_ORIENTATION_VECTOR))
+   VIEW_ANGLE_FACTOR = 0._EB
    DO IO=1,N_ORIENTATION_VECTOR
+      DLO = 0._EB
       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)
+         COSINE_ARRAY(N) = ORIENTATION_VECTOR(1,IO)*DLANG(1,N) + &
+                           ORIENTATION_VECTOR(2,IO)*DLANG(2,N) + &
+                           ORIENTATION_VECTOR(3,IO)*DLANG(3,N)
+         IF (-COSINE_ARRAY(N) > COS_HALF_VIEW_ANGLE(IO)) &
+            DLO = DLO - (ORIENTATION_VECTOR(1,IO)*DLX(N) + ORIENTATION_VECTOR(2,IO)*DLY(N) + ORIENTATION_VECTOR(3,IO)*DLZ(N))
       ENDDO
       NEAREST_RADIATION_ANGLE(IO) = MINLOC(COSINE_ARRAY,DIM=1)
-      VIEW_ANGLE_AREA(IO) = PI/VIEW_ANGLE_AREA(IO)
+      VIEW_ANGLE_FACTOR(IO) = PI/DLO
    ENDDO
    DEALLOCATE(COSINE_ARRAY)
 ENDIF
@@ -3465,7 +3467,7 @@ SUBROUTINE RADIATION_FVM
 USE PHYSICAL_FUNCTIONS, ONLY : GET_VOLUME_FRACTION, GET_MASS_FRACTION
 REAL(EB) :: RAP, AX, AXU, AXD, AY, AYU, AYD, AZ, AZU, AZD, VC, RU, RD, RP, AFD, &
             ILXU, ILYU, ILZU, QVAL, BBF, BBFA, NCSDROP, RSA_RAT,EFLUX,SOOT_MASS_FRACTION, &
-            AIU_SUM,A_SUM,VOL,VC1,AY1,AZ1,COS_DL,AILFU, &
+            AIU_SUM,A_SUM,VOL,VC1,AY1,AZ1,DLO,COS_DLO,AILFU, &
             RAD_Q_SUM_PARTIAL,KFST4_SUM_PARTIAL,ALPHA_CC
 
 INTEGER  :: N,NN,IIG,JJG,KKG,I,J,K,IW,ICF,II,JJ,KK,IOR,IC,IWUP,IWDOWN, &
@@ -4404,17 +4406,17 @@ SUBROUTINE RADIATION_FVM
                                            +TWO_EPSILON_EB)
                      ENDIF
                   ENDIF
-                  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))
+                  COS_DLO = -DOT_PRODUCT(TEMP_ORIENTATION(1:3),DLANG(1:3,N))
+                  IF (COS_DLO > COS_HALF_VIEW_ANGLE(LP%ORIENTATION_INDEX)) THEN
+                     DLO = -(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) * VIEW_ANGLE_AREA(LP%ORIENTATION_INDEX)
+                        BR%BAND(IBND)%ILW(N) = DLO * IL(BC%IIG,BC%JJG,BC%KKG) * VIEW_ANGLE_FACTOR(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) * VIEW_ANGLE_AREA(LP%ORIENTATION_INDEX)
+                        BR%BAND(IBND)%ILW(N) = DLO * IL_UP(BC%IIG,BC%JJG,BC%KKG) * VIEW_ANGLE_FACTOR(LP%ORIENTATION_INDEX)
                      ENDIF
                   ENDIF
                ENDDO PARTICLE_RADIATION_LOOP

Source/read.f90

@@ -97,9 +97,9 @@ SUBROUTINE READ_DATA(DT)
 
 N_ORIENTATION_VECTOR = 0
 ALLOCATE(ORIENTATION_VECTOR(3,0:10))
-ALLOCATE(ORIENTATION_VIEW_ANGLE(0:10))
+ALLOCATE(COS_HALF_VIEW_ANGLE(0:10))
 ORIENTATION_VECTOR(1:3,0) = (/0._EB,0._EB,-1._EB/)
-ORIENTATION_VIEW_ANGLE = 0._EB
+COS_HALF_VIEW_ANGLE = 0._EB
 
 ! Open the input file
 
@@ -5973,10 +5973,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)
+         COS_HALF_VIEW_ANGLE => REALLOCATE(COS_HALF_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
+      COS_HALF_VIEW_ANGLE(N_ORIENTATION_VECTOR) = 0._EB
    ENDIF
    LPC%FREE_AREA_FRACTION = FREE_AREA_FRACTION
    LPC%POROUS_VOLUME_FRACTION = POROUS_VOLUME_FRACTION
@@ -13330,7 +13330,7 @@ SUBROUTINE READ_DEVC
          IF (ABS(ORIENTATION(1)-ORIENTATION_VECTOR(1,I))<TWO_EPSILON_EB .AND. &
              ABS(ORIENTATION(2)-ORIENTATION_VECTOR(2,I))<TWO_EPSILON_EB .AND. &
              ABS(ORIENTATION(3)-ORIENTATION_VECTOR(3,I))<TWO_EPSILON_EB .AND. &
-             ABS(ORIENTATION_VIEW_ANGLE(I)-VIEW_ANGLE)<TWO_EPSILON_EB) THEN
+             ABS(COS_HALF_VIEW_ANGLE(I)-VIEW_ANGLE)<TWO_EPSILON_EB) THEN
             NEW_ORIENTATION_VECTOR = .FALSE.
             ORIENTATION_INDEX = I
             EXIT
@@ -13344,14 +13344,14 @@ SUBROUTINE READ_DEVC
       N_ORIENTATION_VECTOR = N_ORIENTATION_VECTOR + 1
       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)
+         COS_HALF_VIEW_ANGLE => REALLOCATE(COS_HALF_VIEW_ANGLE,0,N_ORIENTATION_VECTOR+10)
       ENDIF
       IF (ALL(ABS(ORIENTATION(1:3))<TWO_EPSILON_EB)) THEN
          WRITE(MESSAGE,'(3A)') 'ERROR(884): DEVC ',TRIM(ID),' components of ORIENTATION are all zero.'
          CALL SHUTDOWN(MESSAGE) ; RETURN
       ENDIF
       ORIENTATION_VECTOR(1:3,N_ORIENTATION_VECTOR) = ORIENTATION(1:3) / NORM2(ORIENTATION)
-      ORIENTATION_VIEW_ANGLE(N_ORIENTATION_VECTOR) = 0._EB
+      COS_HALF_VIEW_ANGLE(N_ORIENTATION_VECTOR) = 0._EB
       ORIENTATION_INDEX = N_ORIENTATION_VECTOR
    ENDIF
 
@@ -14463,7 +14463,7 @@ SUBROUTINE PROC_DEVC
 
       IF (DV%PROP_INDEX > 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)
+            COS_HALF_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