Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer, | intent(in) | :: | IWin | |||
integer, | intent(in) | :: | WinEl | |||
integer, | intent(in) | :: | IHR | |||
integer, | intent(in) | :: | ZoneNum | |||
integer, | intent(in) | :: | iRefPoint | |||
integer, | intent(in) | :: | CalledFrom | |||
integer, | intent(in), | optional | :: | MapNum |
Nodes of different colours represent the following:
Solid arrows point from a procedure to one which it calls. Dashed arrows point from an interface to procedures which implement that interface. This could include the module procedures in a generic interface or the implementation in a submodule of an interface in a parent module. Where possible, edges connecting nodes are given different colours to make them easier to distinguish in large graphs.
Nodes of different colours represent the following:
Solid arrows point from a procedure to one which it calls. Dashed arrows point from an interface to procedures which implement that interface. This could include the module procedures in a generic interface or the implementation in a submodule of an interface in a parent module. Where possible, edges connecting nodes are given different colours to make them easier to distinguish in large graphs.
SUBROUTINE DayltgInterReflectedIllumComplexFenestration(IWin, WinEl, IHR, ZoneNum, iRefPoint, CalledFrom, MapNum)
! SUBROUTINE INFORMATION:
! AUTHOR Simon Vidanovic
! DATE WRITTEN April 2013
! MODIFIED na
! RE-ENGINEERED na
! PURPOSE OF THIS SUBROUTINE:
! Called from CalcDayltgCoefficients for each complex (bsdf) fenestration and reference point in a daylit
! space, for each sun position. Calculates illuminance (EINTSK and EINTSU) at reference point due
! to internally reflected light by integrating to determine the amount of flux from
! sky and ground (and beam reflected from obstructions) transmitted through
! the center of the window and then reflecting this
! light from the inside surfaces of the space.
! METHODOLOGY EMPLOYED: na
! REFERENCES:
!
! USE STATEMENTS:
implicit none ! Enforce explicit typing of all variables in this routine
! SUBROUTINE ARGUMENT DEFINITIONS:
integer, intent(in) :: IWin ! Window index
integer, intent(in) :: WinEl ! Current window element counter
integer, intent(in) :: IHR ! Hour of day
integer, intent(in) :: ZoneNum ! Zone number
integer, intent(in) :: iRefPoint ! reference point counter
integer, intent(in) :: CalledFrom
integer, intent(in), optional :: MapNum
real(r64), dimension(:,:), allocatable :: FLSK ! Sky related luminous flux
real(r64), dimension(:), allocatable :: FLSU ! Sun related luminous flux, excluding entering beam
real(r64), dimension(:), allocatable :: FLSUdisk ! Sun related luminous flux, due to entering beam
real(r64), dimension(:,:), allocatable :: FirstFluxSK ! Sky related first reflected flux
real(r64), dimension(:), allocatable :: FirstFluxSU ! Sun related first reflected flux, excluding entering beam
real(r64), dimension(:), allocatable :: FirstFluxSUdisk ! Sun related first reflected flux, due to entering beam
real(r64), dimension(:,:), allocatable :: ElementLuminanceSky ! sky related luminance at window element (exterior side)
real(r64), dimension(:), allocatable :: ElementLuminanceSun ! sun related luminance at window element (exterior side),
! exluding beam
real(r64), dimension(:), allocatable :: ElementLuminanceSunDisk ! sun related luminance at window element (exterior side),
! due to sun beam
! Total transmitted flux
real(r64) :: FLSKTot(4)
real(r64) :: FLSUTot
real(r64) :: FLSUdiskTot
! Total for first relflected fluxes
real(r64) :: FFSKTot(4)
real(r64) :: FFSUTot
real(r64) :: FFSUdiskTot
real(r64) :: COSIncSun ! cosine of sun incidence angle (from basis elements)
integer :: iSky ! Sky type index: 1=clear, 2=clear turbid, 3=intermediate, 4=overcast
integer :: iConst ! Construction number
integer :: CurCplxFenState
integer :: NIncBasis
integer :: NTrnBasis
integer :: SolBmIndex ! index of current sun position
integer :: iIncElem ! incoming direction counter
integer :: iBackElem ! outgoing direction counter
real(r64) :: LambdaInc ! current lambda value for incoming direction
!real(r64) :: LambdaTrn ! current lambda value for incoming direction
real(r64) :: dirTrans ! directional bsdf transmittance
real(r64) :: ZoneInsideSurfArea
real(r64) :: TransmittedFlux
CurCplxFenState = SurfaceWindow(IWin)%ComplexFen%CurrentState
iConst = SurfaceWindow(IWin)%ComplexFen%State(CurCplxFenState)%Konst
NTrnBasis = ComplexWind(IWin)%Geom(CurCplxFenState)%Trn%NBasis
if (.not.allocated(FLSK)) allocate(FLSK(NTrnBasis, 4))
FLSK = 0.0d0
if (.not.allocated(FLSU)) allocate(FLSU(NTrnBasis))
FLSU = 0.0d0
if (.not.allocated(FLSUdisk)) allocate(FLSUdisk(NTrnBasis))
FLSUdisk = 0.0d0
if (.not.allocated(FirstFluxSK)) allocate(FirstFluxSK(NTrnBasis, 4))
FirstFluxSK = 0.0d0
if (.not.allocated(FirstFluxSU)) allocate(FirstFluxSU(NTrnBasis))
FirstFluxSU = 0.0d0
if (.not.allocated(FirstFluxSUdisk)) allocate(FirstFluxSUdisk(NTrnBasis))
FirstFluxSUdisk = 0.0d0
NIncBasis = ComplexWind(IWin)%Geom(CurCplxFenState)%Inc%NBasis
if (.not.allocated(ElementLuminanceSky)) allocate(ElementLuminanceSky(NIncBasis, 4))
ElementLuminanceSky = 0.0d0
if (.not.allocated(ElementLuminanceSun)) allocate(ElementLuminanceSun(NIncBasis))
ElementLuminanceSun = 0.0d0
if (.not.allocated(ElementLuminanceSunDisk)) allocate(ElementLuminanceSunDisk(NIncBasis))
ElementLuminanceSunDisk = 0.0d0
! Integration over sky/ground/sun elements is done over window incoming basis element and flux is calculated for each
! outgoing direction. This is used to calculate first reflected flux
call ComplexFenestrationLuminances(IWIN, WinEl, NIncBasis, IHR, iRefPoint, ElementLuminanceSky, ElementLuminanceSun, &
ElementLuminanceSunDisk, CalledFrom, MapNum)
! luminance from sun disk needs to include fraction of sunlit area
SolBmIndex = ComplexWind(IWin)%Geom(CurCplxFenState)%SolBmIndex(IHR, timestep)
if (SolBmIndex > 0) then
COSIncSun = ComplexWind(IWin)%Geom(CurCplxFenState)%CosInc(SolBmIndex)
else
COSIncSun = 0.0d0
end if
ElementLuminanceSunDisk = ElementLuminanceSunDisk * SunLitFracHR(IWin, IHR) * COSIncSun
FLSKTot = 0.0d0
FLSUTot = 0.0d0
FLSUdiskTot = 0.0d0
FFSKTot = 0.0d0
FFSUTot = 0.0d0
FFSUdiskTot = 0.0d0
! now calculate flux into each outgoing direction by integrating over all incoming directions
do iBackElem = 1, NTrnBasis
do iIncElem = 1, NIncBasis
LambdaInc = ComplexWind(IWin)%Geom(CurCplxFenState)%Inc%Lamda(iIncElem)
dirTrans = Construct(iConst)%BSDFInput%VisFrtTrans(iIncElem, iBackElem)
do iSky = 1, 4
FLSK(iBackElem, iSky) = FLSK(iBackElem, iSky) + dirTrans * LambdaInc * ElementLuminanceSky(iIncElem, iSky)
end do
FLSU(iBackElem) = FLSU(iBackElem) + dirTrans * LambdaInc * ElementLuminanceSun(iIncElem)
FLSUdisk(iBackElem) = FLSUdisk(iBackElem) + dirTrans * LambdaInc * ElementLuminanceSunDisk(iIncElem)
end do
do iSky = 1, 4
FirstFluxSK(iBackElem, iSky) = FLSK(iBackElem, iSky) * ComplexWind(IWin)%Geom(CurCplxFenState)%AveRhoVisOverlap(iBackElem)
FFSKTot(iSky) = FFSKTot(iSky) + FirstFluxSK(iBackElem, iSky)
FLSKTot(iSky) = FLSKTot(iSky) + FLSK(iBackElem, iSky)
end do
FirstFluxSU(iBackElem) = FLSU(IBackElem) * ComplexWind(IWin)%Geom(CurCplxFenState)%AveRhoVisOverlap(iBackElem)
FFSUTot = FFSUTot + FirstFluxSU(iBackElem)
FLSUTot = FLSUTot + FLSU(iBackElem)
FirstFluxSUdisk(iBackElem) = FLSUdisk(IBackElem) * ComplexWind(IWin)%Geom(CurCplxFenState)%AveRhoVisOverlap(iBackElem)
FFSUdiskTot = FFSUdiskTot + FirstFluxSUdisk(iBackElem)
FLSUDiskTot = FLSUDiskTot + FLSUDisk(iBackElem)
end do
ZoneInsideSurfArea = ZoneDaylight(ZoneNum)%TotInsSurfArea
do iSky = 1, 4
EINTSK(iSky, 1, IHR) = FFSKTot(iSky) * (Surface(IWin)%Area/SurfaceWindow(IWin)%GlazedFrac) &
/ (ZoneInsideSurfArea * (1.0d0 - ZoneDaylight(ZoneNum)%AveVisDiffReflect))
end do
EINTSU(1,IHR) = FFSUTot * (Surface(IWin)%Area/SurfaceWindow(IWin)%GlazedFrac) &
/ (ZoneInsideSurfArea * (1.0d0 - ZoneDaylight(ZoneNum)%AveVisDiffReflect))
EINTSUdisk(1,IHR) = FFSUdiskTot * (Surface(IWin)%Area/SurfaceWindow(IWin)%GlazedFrac) &
/ (ZoneInsideSurfArea * (1.0d0 - ZoneDaylight(ZoneNum)%AveVisDiffReflect))
if (allocated(FLSK)) deallocate(FLSK)
if (allocated(FLSU)) deallocate(FLSU)
if (allocated(FLSUdisk)) deallocate(FLSUdisk)
if (allocated(FirstFluxSK)) deallocate(FirstFluxSK)
if (allocated(FirstFluxSU)) deallocate(FirstFluxSU)
if (allocated(FirstFluxSUdisk)) deallocate(FirstFluxSUdisk)
if (allocated(ElementLuminanceSky)) deallocate(ElementLuminanceSky)
if (allocated(ElementLuminanceSun)) deallocate(ElementLuminanceSun)
if (allocated(ElementLuminanceSunDisk)) deallocate(ElementLuminanceSunDisk)
END SUBROUTINE DayltgInterReflectedIllumComplexFenestration