SUBROUTINE CalcComplexWindowThermal(SurfNum,ConstrNum,HextConvCoeff,SurfInsideTemp,SurfOutsideTemp,SurfOutsideEmiss,CalcCondition)
! SUBROUTINE INFORMATION:
! AUTHOR B. Griffith
! DATE WRITTEN October 2009
! MODIFIED Simon Vidanovic
! RE-ENGINEERED September 2011
! PURPOSE OF THIS SUBROUTINE:
! wrapper between E+ and TARCOG
! METHODOLOGY EMPLOYED:
! draft out an attempt for proof-of-concept, to reuse native TARCOG implementation
! based off of 1-26-2009 version of WinCOG/TARCOG solution from Carli, Inc.
! REFERENCES:
! na
! USE STATEMENTS:
!USE WindowTARCOGManager, ONLY: tarcog
USE DataBSDFWindow
USE DataZoneEquipment, ONLY : ZoneEquipConfig
USE DataLoopNode, ONLY : Node
USE Psychrometrics, ONLY : PsyCpAirFnWTdb,PsyTdpFnWPb
USE General, ONLY : InterpSlatAng , InterpSw ! Function for slat angle interpolation
USE InputProcessor, ONLY : SameString
USE DataHeatBalSurface, ONLY : HcExtSurf
USE DataGlobals, ONLY: StefanBoltzmann
USE TARCOGGassesParams, ONLY : maxgas
USE TARCOGParams, ONLY : maxlay, maxlay1
USE DataHeatBalance, ONLY : GasCoeffsAir, SupportPillar
USE TARCOGMain, ONLY : TARCOG90
IMPLICIT NONE ! Enforce explicit typing of all variables in this routine
! SUBROUTINE ARGUMENT DEFINITIONS:
INTEGER, INTENT(IN) :: SurfNum ! Surface number
INTEGER, INTENT(IN) :: CalcCondition ! Calucation condition (summer, winter or no condition)
INTEGER, INTENT(INOUT) :: ConstrNum ! Construction number
REAL(r64), INTENT(IN) :: HextConvCoeff ! Outside air film conductance coefficient
REAL(r64), INTENT(INOUT) :: SurfInsideTemp ! Inside window surface temperature
REAL(r64), INTENT(INOUT) :: SurfOutsideTemp ! Outside surface temperature (C)
! (temperature of innermost face) [C]
REAL(r64), INTENT(INOUT) :: SurfOutsideEmiss
!INTEGER, INTENT(IN) :: CurrentThermalModelNumber
! SUBROUTINE PARAMETER DEFINITIONS:
!INTEGER, PARAMETER :: maxlay = 100 ! maximum number of layers (including laminates)
!INTEGER, PARAMETER :: maxgas = 10 ! maximum number of individual gasses
!INTEGER, PARAMETER :: maxlay1 = maxlay+1 ! maximum number of 'gaps', including in and out (maxlay+1)
!REAL(r64), PARAMETER :: StefanBoltzmannConst = 5.6697d-8 ! Stefan-Boltzmann constant in W/(m2*K4)
! INTERFACE BLOCK SPECIFICATIONS:
! na
! DERIVED TYPE DEFINITIONS:
! na
! SUBROUTINE LOCAL VARIABLE DECLARATIONS:
! TARCOG Inputs:
INTEGER :: nlayer = 0 ! Number of glazing layers
INTEGER :: iwd = 0 ! Wind direction: 0 - windward, 1 - leeward
REAL(r64) :: tout = 0.0D0 ! Outdoor temperature [K]
REAL(r64) :: tind = 0.0D0 ! Indoor temperature [K]
REAL(r64) :: trmin = 0.0D0 ! Indoor mean radiant temperature [K]
REAL(r64) :: wso = 0.0D0 ! Outdoor wind speed [m/s]
REAL(r64) :: wsi = 0.0D0 ! Inside forced air speed [m/s]
REAL(r64) :: dir = 0.0D0 ! Direct solar radiation [W/m^2]
INTEGER :: isky = 0 ! Flag for sky temperature (Tsky) and sky emittance (esky)
! 0 - both tsky and esky are specified
! 1 - tsky specified, esky = 1
! 2 - Swinbank model for effective sky emittance
REAL(r64) :: tsky = 0.0D0 ! Night sky temperature [K]
REAL(r64) :: esky = 0.0D0 ! Effective night sky emittance
REAL(r64) :: fclr = 0.0D0 ! Fraction of sky that is clear
REAL(r64) :: VacuumPressure ! maximal pressure for gas to be considered as vacuum [Pa]
REAL(r64) :: VacuumMaxGapThickness ! maximal gap thickness for which vacuum calculation will work without issuing
! warning message
REAL(r64) , DIMENSION(maxlay) :: gap = 0.0D0 ! Vector of gap widths [m] {maxlay}
REAL(r64) , DIMENSION(maxlay) :: thick = 0.0D0 ! Vector of glass thicknesses [m] {maxlay}
REAL(r64) , DIMENSION(maxlay) :: scon = 0.0D0 ! Vector of conductivities of each glazing layer [W/m.K] {maxlay}
REAL(r64) , DIMENSION(maxlay * 2 ) :: tir = 0.0D0 ! Vector of IR transmittances of each layer {2*maxlay - 2 surfaces per layer}
REAL(r64) , DIMENSION(maxlay * 2 ) :: emis = 0.0D0 ! Vector of IR emittances of each surface {2*maxlay - 2 surfaces per layer}
INTEGER, DIMENSION(maxlay) :: SupportPlr = 0 ! Shows whether or not gap have support pillar
! 0 - does not have support pillar
! 1 - have support pillar
REAL(r64), DIMENSION(maxlay) :: PillarSpacing = 0.0d0 ! Pillar spacing for each gap (used in case there is support pillar)
REAL(r64), DIMENSION(maxlay) :: PillarRadius = 0.0d0 ! Pillar radius for each gap (used in case there is support pillar)
REAL(r64) :: totsol = 0.0D0 ! Total solar transmittance of the IGU
REAL(r64) :: tilt = 0.0D0 ! Window tilt [degrees]
REAL(r64) , DIMENSION(maxlay) :: asol = 0.0D0 ! Vector of Absorbed solar energy fractions for each layer {maxlay}
REAL(r64) :: height = 0.0D0 ! IGU cavity height [m]
REAL(r64) :: heightt = 0.0D0 ! Total window height [m]
REAL(r64) :: width = 0.0D0 ! Window width [m]
REAL(r64) , DIMENSION(maxlay+1) :: presure = 0.0D0 ! Vector of gas pressures in gaps [N/m^2] {maxlay+1}
!Deflection
!Tarcog requires deflection as input parameters. Deflection is NOT used in EnergyPlus simulations
INTEGER :: CalcDeflection ! Deflection calculation flag:
! 0 - no deflection calculations
! 1 - perform deflection calculation (input is Pressure/Temp)
! 2 - perform deflection calculation (input is measured deflection)
REAL(r64) :: Pa ! Atmospheric (outside/inside) pressure (used onlu if CalcDeflection = 1)
REAL(r64) :: Pini ! Initial presssure at time of fabrication (used only if CalcDeflection = 1)
REAL(r64) :: Tini ! Initial temperature at time of fabrication (used only if CalcDeflection = 1)
REAL(r64), DIMENSION(maxlay-1) :: GapDefMax ! Vector of gap widths in deflected state. It will be used as input
! if CalcDeflection = 2. In case CalcDeflection = 1 it will return recalculated
! gap widths. [m]
REAL(r64), DIMENSION(maxlay) :: YoungsMod ! Vector of Young's modulus. [m]
REAL(r64), DIMENSION(maxlay) :: PoissonsRat ! Vector of Poisson's Ratios. [m]
REAL(r64), DIMENSION(maxlay) :: LayerDef ! Vector of layers deflection. [m]
INTEGER, DIMENSION(maxlay+1, maxgas) :: iprop = 1 ! Matrix of gas codes - see above {maxgap x maxgas}
REAL(r64), DIMENSION(maxlay+1, maxgas) :: frct = 0.0D0 ! Matrix of mass percentages in gap mixtures {maxgap x maxgas}
REAL(r64), DIMENSION(maxgas, 3) :: gcon = 0.0D0 ! Matrix of constants for gas conductivity calc
! (A, B, C for max of 10 gasses) {maxgas x 3}
REAL(r64), DIMENSION(maxgas, 3) :: gvis = 0.0D0 ! Matrix of constants for gas dynamic viscosity calc
! (A, B, C for max of 10 gasses) {maxgas x 3}
REAL(r64), DIMENSION(maxgas, 3) :: gcp = 0.0D0 ! Matrix of constants for gas specific heat calc at constant pressure
! (A, B, C for max of 10 gasses) {maxgas x 3}
REAL(r64), DIMENSION(maxgas) :: wght = 0.0D0 ! Vector of Molecular weights for gasses {maxgas}
REAL(r64), DIMENSION(maxgas) :: gama = 0.0D0 ! Vector of spefic heat ration for low pressure calc {maxgas}
LOGICAL :: feedData = .FALSE. !flag to notify if data needs to be feed into gas arrays
INTEGER, DIMENSION(maxlay + 1) :: nmix = 0 ! Vector of number of gasses in gas mixture of each gap {maxlay+1}
REAL(r64) :: hin = 0.0D0 ! Indoor combined film coefficient (if non-zero) [W/m^2.K]
REAL(r64) :: hout = 0.0D0 ! Outdoor combined film coefficient (if non-zero) [W/m^2.K]
INTEGER, DIMENSION(2) :: ibc = 0 ! Vector of boundary condition flags (ibc(1) - outdoor, ibc(2) - indoor)
! 0 - h to be calculated;
! 1 - combined film coefficient (h) prescribed;
! 2 - convective film coefficient (hc) prescribed.
! Also used in old algorithms for calculating h, accessible through
! negative values for flags:
! -1 - old SPC142 correlation
! -2 - Klems-Yazdanian correlation (applicable to outdoor only)
! -3 - Kimura correlation (applicable to outdoor only)
REAL(r64), DIMENSION(maxlay) :: Atop = 0.0D0 ! Vector with areas of top openings - between SD layers and top of
! glazing cavity, for each layer [m^2] {maxlay} *
REAL(r64), DIMENSION(maxlay) :: Abot = 0.0D0 ! Vector with areas of bottom openings - between SD layers
! and bottom of glazing cavity [m^2] {maxlay}
REAL(r64), DIMENSION(maxlay) :: Al = 0.0D0 ! Vector with areas of left-hand side openings - between SD layers
! and left end of glazing cavity [m^2] {maxlay}
REAL(r64), DIMENSION(maxlay) :: Ar = 0.0D0 ! Vector of areas of right-hand side openings - between SD layers
! and right end of glazing cavity [m^2] {maxlay}
REAL(r64), DIMENSION(maxlay) :: Ah = 0.0D0 ! Vector of total areas of holes for each SD [m^2] {maxlay}
REAL(r64), DIMENSION(maxlay) :: SlatThick = 0.0D0 ! Thickness of the slat material [m] {maxlay} **
REAL(r64), DIMENSION(maxlay) :: SlatWidth = 0.0D0 ! Slat width [m] {maxlay}
REAL(r64), DIMENSION(maxlay) :: SlatAngle = 0.0D0 ! Slat tilt angle [deg] {maxlay}
REAL(r64), DIMENSION(maxlay) :: SlatCond = 0.0D0 ! Conductivity of the slat material [W/m.K] {maxlay}
REAL(r64), DIMENSION(maxlay) :: SlatSpacing = 0.0D0 ! Distance between slats [m] {maxlay}
REAL(r64), DIMENSION(maxlay) :: SlatCurve = 0.0D0 ! Curvature radius of the slat [m] {maxlay}
REAL(r64), DIMENSION(maxlay+1) :: vvent = 0.0D0 ! Vector of velocities for forced ventilation, for each gap, and for
! outdoor and indoor environment [m/s] {maxlay+1} ***
REAL(r64), DIMENSION(maxlay+1) :: tvent = 0.0D0 ! Vector of temperatures of ventilation gas for forced ventilation, for each
! gap, and for outdoor and indoor environment [K] {maxlay+1}
INTEGER, DIMENSION(maxlay) :: LayerType = 0 ! Glazing layer type flag {maxlay}:
! 0 - Specular layer,
! 1 - Venetian blind (SD)
! 2 - Woven shade (SD) (not implemented)
! 3 - Diffuse shade (not implemented)
INTEGER, DIMENSION(maxlay) :: nslice = 0 ! Vector of numbers of slices in a laminated glazing layers
! (0 - monolithic layer) {maxlay}
REAL(r64), DIMENSION(maxlay) :: LaminateA = 0.0D0 ! Left-hand side array for creating slice equations {maxlay}
REAL(r64), DIMENSION(maxlay) :: LaminateB = 0.0D0 ! Right-hand side array for creating slice equations {maxlay}
REAL(r64), DIMENSION(maxlay) :: sumsol = 0.0D0 ! Array of absorbed solar energy fractions for each laminated
! glazing layer [W/m^2] {maxlay}
INTEGER :: standard = 1 ! Calculation standard switch:
! 1 - ISO 15099,
! 2 - EN673 / ISO 10292 Declared,
! 3 - EN673 / ISO 10292 Design.
INTEGER :: ThermalMod = 0 ! Thermal model:
! 0 - ISO15099
! 1 - Scaled Cavity Width (SCW)
! 2 - Convective Scalar Model (CSM)
INTEGER :: Debug_mode = 0 ! Switch for debug output files:
! 0 - donÂ’'t create debug output files
! 1 - append results to existing debug output file (where applicable)
! 2 - store results in a new debug output file
CHARACTER(len=256) :: Debug_dir = ' ' ! Target directory for debug files (pointer to a character array)
CHARACTER(len=256) :: Debug_file = 'Test' ! Template file name used to create debug output files
INTEGER(4) :: Window_ID = -1 ! ID of the window (long integer value, passed by W6)
INTEGER(4) :: IGU_ID = -1 ! ID of the IGU (long integer value, passed by W6)
REAL(r64) :: SDScalar = 0.0D0 ! SD convection factor (value between 0 and 1)
! 0.0 - No SD layer
! 1.0 - Closed SD
!
! Notes: * vvent, tvent, Atop, Abot, Al, Ar and Ah are considered for SD layers only.
! ** SlatThick, SlatWidth, SlatAngle, SlatCond, SlatSpacing, SlatCurve
! are used for Venetian blind layers only.
! *** For vvent & tvent: vvent(1) - exterior, vvent(nlayer+1) - interior.
! **** Forced ventilation calculation is not active at this time.
!
! TARCOG Output:
REAL(r64) , DIMENSION(maxlay * 2 ) :: theta = 0.0D0 ! Vector of average temperatures of glazing surfaces [K] {2*maxlay}
REAL(r64) , DIMENSION(maxlay * 2 + 1) :: q = 0.0D0 ! Vector of various heat fluxes [W/m^2] {2*maxlay+1},
! depending on element index:
! 1 = qout (heat flux from outer-most glazing surface to outdoor space)
! 2*i = qpane(i) (heat flux through i-th glazing layer)
! 2*i-1 = qgap(i) (heat flux from i-th glazing cavity to indoor-faced
! surface of the adjacent glazing layer)
! 2*nlayer+1 = qin (heat flux from indoor space to inner-most glazing
! surface)
REAL(r64) , DIMENSION(maxlay1) :: qprim = 0.0D0 ! Vector of heat fluxes from the outdoor-faced surfaces of glazing layers
! towards the adjacent glazing cavity [W/m2]
REAL(r64) , DIMENSION(maxlay1) :: qv = 0.0D0 ! Vector of heat fluxes to each gap by ventillation [W/m^2]
REAL(r64) :: ufactor = 0.0D0 ! Center of glass U-value [W/m^2.K]
REAL(r64) :: sc = 0.0D0 ! Shading Coefficient
REAL(r64) :: hflux = 0.0D0 ! Net heat flux between room and window [W/m^2]
REAL(r64) :: rhum = 0.0D0 ! Relative humidity at inner-most surface for condensation
REAL(r64) :: rhout = 0.0D0 ! Relative humidity at outer-most surface for condensation
REAL(r64) :: hcin = 0.0D0 ! Indoor convective surface heat transfer coefficient [W/m^2.K]
REAL(r64) :: hcout = 0.0D0 ! Outdoor convective surface heat transfer coefficient [W/m^2.K]
REAL(r64) :: hrin = 0.0D0 ! Indoor radiative surface heat transfer coefficient [W/m^2.K]
REAL(r64) :: hrout = 0.0D0 ! Outdoor radiative surface heat transfer coefficient [W/m^2.K]
REAL(r64) , DIMENSION(maxlay1) :: hcgap = 0.0D0 ! Convective part of gap effective conductivity {maxlay}
REAL(r64) , DIMENSION(maxlay1) :: hrgap = 0.0D0 ! Radiative part of gap effective conductivity (including in and out)
REAL(r64) :: shgc = 0.0D0 ! Solar heat gain coefficient - per ISO 15099
REAL(r64) :: shgct = 0.0D0 ! Solar heat gain coefficient - per old procedure
REAL(r64) :: tamb = 0.0D0 ! Outdoor environmental temperature [K]
REAL(r64) :: troom = 0.0D0 ! Indoor environmental temperature [K]
REAL(r64) , DIMENSION(maxlay) :: hg = 0.0D0 ! Gas conductance of the glazing cavity
! [W/m^2.K] - EN673 and ISO 10292 procedure
REAL(r64) , DIMENSION(maxlay) :: hr = 0.0D0 ! Radiation conductance of the glazing cavity
! [W/m^2.K] - EN673 and ISO 10292 procedure
REAL(r64) , DIMENSION(maxlay) :: hs = 0.0D0 ! Thermal conductance of the glazing cavity
! [W/m^2.K] - EN673 and ISO 10292 procedure
REAL(r64) :: he = 0.0D0 ! External heat transfer coefficient [W/m^2.K] - EN673 and ISO 10292 procedure
REAL(r64) :: hi = 0.0D0 ! Internal heat transfer coefficient [W/m^2.K] - EN673 and ISO 10292 procedure
REAL(r64) , DIMENSION(maxlay + 1) :: Ra = 0.0D0 ! Vector of Rayleigh numbers, for each gap {maxlay}
REAL(r64) , DIMENSION(maxlay + 1) :: Nu = 0.0D0 ! Vector of Nusselt numbers, for each gap {maxlay}
INTEGER :: nperr = 0 ! Error code
REAL(r64) :: ShadeEmisRatioOut = 0.0D0 ! Ratio of modified to glass emissivity at the outermost glazing surface
REAL(r64) :: ShadeEmisRatioIn = 0.0D0 ! Ratio of modified to glass emissivity at the innermost glazing surface
REAL(r64) :: ShadeHcRatioOut = 0.0D0 ! Ratio of modified to unshaded Hc at the outermost glazing surface
REAL(r64) :: ShadeHcRatioIn = 0.0D0 ! Ratio of modified to unshaded Hc at the innermost glazing surface
REAL(r64) :: HcUnshadedOut = 0.0D0 ! Hc value at outdoor surface of an unshaded subsystem [W/m^2.K]
REAL(r64) :: HcUnshadedIn = 0.0D0 ! Hc value at indoor surface of an unshaded subsystem [W/m^2.K]
REAL(r64) , DIMENSION(maxlay) :: Keff = 0.0D0 ! Vector of keff values for gaps [W/m.K] {maxlay}
REAL(r64) , DIMENSION(maxlay-1) :: ShadeGapKeffConv = 0.0D0 ! Vector of convective keff values for areas above/below
! SD layers [W/m.K] {maxlay-1}
INTEGER :: ZoneNum ! Zone number corresponding to SurfNum
REAL(r64) :: locTCSpecTemp ! The temperature corresponding to the specified optical properties of the TC layer
REAL(r64) :: locTCLayerTemp ! TC layer temperature at each time step. C
LOGICAL :: locTCFlag =.False. ! True if this surface is a TC window
REAL(r64) :: deltaTemp(100) = 0.0d0
INTEGER :: i
INTEGER :: iMinDT(1) = 0
INTEGER :: IDConst(100) = 0
REAL(r64) :: dT0 = 0.0d0
REAL(r64) :: dT1 = 0.0d0
INTEGER :: IConst ! Construction number
INTEGER :: TotLay ! Total number of layers in a construction
! (sum of solid layers and gap layers)
INTEGER :: Lay ! Layer number
INTEGER :: LayPtr ! Material number for a layer
INTEGER :: ShadingLayPtr ! Shading layer pointer for effective temperature calculations
INTEGER :: GlassLayPtr ! Glass layer pointer for effective temperature calculations
REAL(r64) :: EpsGlassIR
REAL(r64) :: RhoGlassIR
REAL(r64) :: TauShadeIR
REAL(r64) :: EpsShadeIR
REAL(r64) :: RhoShadeIR
INTEGER :: IGlass ! glass layer number (1,2,3,...)
INTEGER :: IGap ! Gap layer number (1,2,...)
INTEGER :: ShadeLayPtr ! Material number corresponding to a shade layer
INTEGER :: TotGlassLay ! Total number of glass layers in a construction
INTEGER :: BlNum ! Blind number
INTEGER :: ZoneEquipConfigNum
INTEGER :: NodeNum
REAL(r64) :: SumSysMCp ! Zone sum of air system MassFlowRate*Cp
REAL(r64) :: SumSysMCpT ! Zone sum of air system MassFlowRate*Cp*T
REAL(r64) :: MassFlowRate
REAL(r64) :: NodeTemp
REAL(r64) :: CpAir
REAL(r64) :: RefAirTemp ! reference air temperatures
INTEGER :: tmpGasType
INTEGER :: k ! Layer counter
INTEGER :: SurfNumAdj ! An interzone surface's number in the adjacent zone
INTEGER :: ZoneNumAdj ! An interzone surface's adjacent zone number
INTEGER :: ShadeFlag ! Flag indicating whether shade or blind is on, and shade/blind position
INTEGER :: imix
REAL(r64) :: EffShBlEmiss ! Effective interior shade or blind emissivity
REAL(r64) :: EffGlEmiss ! Effective inside glass emissivity when interior shade or blind
REAL(r64) :: IncidentSolar ! Solar incident on outside of window (W)
REAL(r64) :: ConvHeatFlowNatural ! Convective heat flow from gap between glass and interior shade or blind (W)
REAL(r64) :: ShadeArea ! shade/blind area (m2)
REAL(r64) :: sconsh ! shade/blind conductance (W/m2-K)
REAL(r64) :: CondHeatGainShade ! Conduction through shade/blind, outside to inside (W)
REAL(r64) :: ShGlReflFacIR ! Factor for long-wave inter-reflection between shade/blind and adjacent glass
REAL(r64) :: RhoGlIR1,RhoGlIR2 ! Long-wave reflectance of glass surface facing shade/blind; 1=exterior shade/blind,
! 2=interior shade/blind
REAL(r64) :: RhoShIR1,RhoShIR2 ! Long-wave reflectance of shade/blind surface facing glass; 1=interior shade/blind,
! 2=exterior shade/blind
REAL(r64) :: EpsShIR1,EpsShIR2 ! Long-wave emissivity of shade/blind surface facing glass; 1=interior shade/blind,
! 2=exterior shade/blind
REAL(r64) :: TauShIR ! Long-wave transmittance of isolated shade/blind
REAL(r64) :: NetIRHeatGainShade ! Net IR heat gain to zone from interior shade/blind (W)
REAL(r64) :: NetIRHeatGainGlass ! Net IR heat gain to zone from shade/blind side of glass when interior
! shade/blind is present. Zero if shade/blind has zero IR transmittance (W)
REAL(r64) :: ConvHeatGainFrZoneSideOfShade ! Convective heat gain to zone from side of interior shade facing zone (W)
REAL(r64) :: ConvHeatGainFrZoneSideOfGlass ! Convective heat gain to zone from side of glass facing zone when
! no interior shade/blind is present (W)
REAL(r64) :: CondHeatGainGlass ! Conduction through inner glass layer, outside to inside (W)
REAL(r64) :: TotAirflowGap ! Total volumetric airflow through window gap (m3/s)
REAL(r64) :: TAirflowGapOutlet ! Temperature of air leaving airflow gap between glass panes (K)
REAL(r64) :: TAirflowGapOutletC ! Temperature of air leaving airflow gap between glass panes (C)
REAL(r64) :: ConvHeatFlowForced ! Convective heat flow from forced airflow gap (W)
REAL(r64) :: InletAirHumRat ! Humidity ratio of air from window gap entering fan
REAL(r64) :: ZoneTemp ! Zone air temperature (C)
REAL(r64) :: CpAirOutlet ! Heat capacity of air from window gap (J/kg-K)
REAL(r64) :: CpAirZone ! Heat capacity of zone air (J/kg-K)
REAL(r64) :: ConvHeatGainToZoneAir ! Convective heat gain to zone air from window gap airflow (W)
INTEGER :: ConstrNumSh ! Construction number with shading device
REAL(r64) :: TransDiff ! Diffuse shortwave transmittance
INTEGER :: CalcSHGC = 0 ! SHGC calculations are not necessary for E+ run
INTEGER :: NumOfIterations = 0 !
integer :: GasType !locally used coefficent to point at correct gas type
integer :: ICoeff
CHARACTER (len=2000) :: tarcogErrorMessage ! store error text from tarcog
!Simon: locally used variables
integer :: ngllayer, nglface, nglfacep, tempInt
integer :: PillarPtr, DeflectionPtr, GasPointer, ThermalModelNum
real(r64) :: Rmir, outir, Ebout ! IR radiance of window's interior surround (W/m2)
REAL(r64) :: dominantGapWidth ! store value for dominant gap width. Used for airflow calculations
! fill local vars
CalcDeflection = 0
CalcSHGC = 0
if (CalcCondition == noCondition) then
ConstrNum = Surface(SurfNum)%Construction
SurfNumAdj = Surface(SurfNum)%ExtBoundCond
ShadeFlag = SurfaceWindow(SurfNum)%ShadingFlag
end if
TotGlassLay = Construct(ConstrNum)%TotGlassLayers
ngllayer = Construct(ConstrNum)%TotGlassLayers
nglface = 2*ngllayer
nglfacep = nglface
hrin = 0.0d0
hcin = 0.0d0
hrout = 0.0d0
hcout = 0.0d0
Pa = OutBaroPress
ThermalModelNum = Construct(ConstrNum)%BSDFInput%ThermalModel
standard = WindowThermalModel(ThermalModelNum)%CalculationStandard
ThermalMod = WindowThermalModel(ThermalModelNum)%ThermalModel
CalcDeflection = WindowThermalModel(ThermalModelNum)%DeflectionModel
SDScalar = WindowThermalModel(ThermalModelNum)%SDScalar
VacuumPressure = WindowThermalModel(ThermalModelNum)%VacuumPressureLimit
Tini = WindowThermalModel(ThermalModelNum)%InitialTemperature - KelvinConv
Pini = WindowThermalModel(ThermalModelNum)%InitialPressure
if (CalcCondition == noCondition) then
ZoneNum = Surface(SurfNum)%Zone
end if
nlayer = Construct(ConstrNum)%TotSolidLayers
isky = 3 ! IR radiation is provided from external source
iwd = 0 ! assume windward for now. TODO compare surface normal with wind direction
if (CalcCondition == noCondition) then
! determine reference air temperature for this surface
SELECT CASE (Surface(SurfNum)%TAirRef)
CASE (ZoneMeanAirTemp)
RefAirTemp = MAT(ZoneNum)
CASE (AdjacentAirTemp)
RefAirTemp = TempEffBulkAir(SurfNum)
CASE (ZoneSupplyAirTemp)
! determine ZoneEquipConfigNum for this zone
! ControlledZoneAirFlag = .FALSE.
ZoneEquipConfigNum = ZoneNum
! DO ZoneEquipConfigNum = 1, NumOfControlledZones
! IF (ZoneEquipConfig(ZoneEquipConfigNum)%ActualZoneNum /= ZoneNum) CYCLE
! ControlledZoneAirFlag = .TRUE.
! EXIT
! END DO ! ZoneEquipConfigNum
! check whether this zone is a controlled zone or not
IF (.NOT. Zone(ZoneNum)%IsControlled) THEN
CALL ShowFatalError('Zones must be controlled for Ceiling-Diffuser Convection model. No system serves zone '// &
TRIM(Zone(ZoneNum)%Name))
RETURN
END IF
! determine supply air conditions
SumSysMCp = 0.0d0
SumSysMCpT = 0.0d0
DO NodeNum = 1, ZoneEquipConfig(ZoneEquipConfigNum)%NumInletNodes
NodeTemp = Node(ZoneEquipConfig(ZoneEquipConfigNum)%InletNode(NodeNum))%Temp
MassFlowRate = Node(ZoneEquipConfig(ZoneEquipConfigNum)%InletNode(NodeNum))%MassFlowRate
CpAir = PsyCpAirFnWTdb(ZoneAirHumRat(ZoneNum), NodeTemp,'CalcComplexWindowThermal')
SumSysMCp = SumSysMCp + MassFlowRate * CpAir
SumSysMCpT = SumSysMCpT + MassFlowRate * CpAir * NodeTemp
END DO
! a weighted average of the inlet temperatures.
RefAirTemp = SumSysMCpT/SumSysMCp
CASE DEFAULT
! currently set to mean air temp but should add error warning here
RefAirTemp = MAT(ZoneNum)
END SELECT
tind = RefAirTemp + KelvinConv ! Inside air temperature
! now get "outside" air temperature
IF(SurfNumAdj > 0) THEN ! Interzone window
ZoneNumAdj = Surface(SurfNumAdj)%Zone
! determine reference air temperature for this surface
SELECT CASE (Surface(SurfNumAdj)%TAirRef)
CASE (ZoneMeanAirTemp)
RefAirTemp = MAT(ZoneNumAdj)
CASE (AdjacentAirTemp)
RefAirTemp = TempEffBulkAir(SurfNumAdj)
CASE (ZoneSupplyAirTemp)
! determine ZoneEquipConfigNum for this zone
ZoneEquipConfigNum = ZoneNum
! check whether this zone is a controlled zone or not
IF (.NOT. Zone(ZoneNum)%IsControlled) THEN
CALL ShowFatalError('Zones must be controlled for Ceiling-Diffuser Convection model. No system serves zone '// &
TRIM(Zone(ZoneNum)%Name))
RETURN
END IF
! determine supply air conditions
SumSysMCp = 0.0d0
SumSysMCpT = 0.0d0
DO NodeNum = 1, ZoneEquipConfig(ZoneEquipConfigNum)%NumInletNodes
NodeTemp = Node(ZoneEquipConfig(ZoneEquipConfigNum)%InletNode(NodeNum))%Temp
MassFlowRate = Node(ZoneEquipConfig(ZoneEquipConfigNum)%InletNode(NodeNum))%MassFlowRate
CpAir = PsyCpAirFnWTdb(ZoneAirHumRat(ZoneNumAdj), NodeTemp, 'CalcComplexWindowThermal')
SumSysMCp = SumSysMCp + MassFlowRate * CpAir
SumSysMCpT = SumSysMCpT + MassFlowRate * CpAir * NodeTemp
END DO
! a weighted average of the inlet temperatures.
RefAirTemp = SumSysMCpT/SumSysMCp
CASE DEFAULT
! currently set to mean air temp but should add error warning here
RefAirTemp = MAT(ZoneNumAdj)
END SELECT
Tout = RefAirTemp + KelvinConv ! outside air temperature
tsky = MRT(ZoneNumAdj) + KelvinConv ! TODO this misses IR from sources such as high temp radiant and baseboards
! ! Add long-wave radiation from adjacent zone absorbed by glass layer closest to the adjacent zone.
!
! AbsRadGlassFace(1) = AbsRadGlassFace(1) + QRadThermInAbs(SurfNumAdj)
!
! ! The IR radiance of this window's "exterior" surround is the IR radiance
! ! from surfaces and high-temp radiant sources in the adjacent zone
!
outir = SurfaceWindow(SurfNumAdj)%IRfromParentZone + QHTRadSysSurf(SurfNumAdj) + QHWBaseboardSurf(SurfNumAdj)
ELSE ! Exterior window (ExtBoundCond = 0)
IF(Surface(SurfNum)%ExtWind) THEN ! Window is exposed to wind (and possibly rain)
IF(IsRain) THEN ! Raining: since wind exposed, outside window surface gets wet
tout = Surface(SurfNum)%OutWetBulbTemp + KelvinConv
ELSE ! Dry
tout = Surface(SurfNum)%OutDryBulbTemp + KelvinConv
END IF
ELSE ! Window not exposed to wind
tout = Surface(SurfNum)%OutDryBulbTemp + KelvinConv
END IF
!tsky = SkyTemp + TKelvin
tsky = SkyTempKelvin
Ebout = sigma * tout**4
outir = Surface(SurfNum)%ViewFactorSkyIR * (AirSkyRadSplit(SurfNum)*sigma*tsky**4 + &
(1.0d0-AirSkyRadSplit(SurfNum))*Ebout) + Surface(SurfNum)%ViewFactorGroundIR * Ebout
END IF
hin = HConvIn(SurfNum) ! Room-side surface convective film conductance
ibc(2) = 0 ! convective coefficient on indoor side will be recalculated (like in Winkelmann routines)
!hcout=HExtConvCoeff ! Exterior convection coefficient is passed in from outer routine
hout=HExtConvCoeff ! Exterior convection coefficient is passed in from outer routine
ibc(1) = 2 ! prescribed convective film coeff on outdoor side
tilt = Surface(SurfNum)%Tilt
height = Surface(SurfNum)%Height
heightt = height ! for now put same window and glazing pocket hights
width = Surface(SurfNum)%Width
!indoor mean radiant temperature.
! IR incident on window from zone surfaces and high-temp radiant sources
rmir = SurfaceWindow(SurfNum)%IRfromParentZone + QHTRadSysSurf(SurfNum) + QHWBaseboardSurf(SurfNum) + &
QSteamBaseboardSurf(SurfNum) + QElecBaseboardSurf(SurfNum)
trmin = ( rmir / StefanBoltzmann)**0.25d0 ! TODO check model equation.
! outdoor wind speed
IF (.NOT. Surface(SurfNum)%ExtWind) THEN
wso = 0.0d0 ! No wind exposure
!ELSE IF (Surface(SurfNum)%Class == SurfaceClass_Window .AND. SurfaceWindow(SurfNum)%ShadingFlag == ExtShadeOn) THEN
! wso = 0.0 ! Assume zero wind speed at outside glass surface of window with exterior shade
ELSE
wso = Surface(SurfNum)%WindSpeed
ENDIF
! indoor wind speed
wsi = 0.0d0 ! assumuption (TODO, what to use for inside air velocity?)
fclr = 1.0D0 - CloudFraction
end if
!now fill layer data
!IConst = ConstrNum
!IF(ShadeFlag==IntShadeOn.OR.ShadeFlag==ExtShadeOn.OR.ShadeFlag==IntBlindOn.OR.ShadeFlag==ExtBlindOn &
! .OR.ShadeFlag==BGShadeOn.OR.ShadeFlag==BGBlindOn.OR.ShadeFlag==ExtScreenOn) THEN
! IConst = Surface(SurfNum)%ShadedConstruction
! IF(Surfacewindow(SurfNum)%StormWinFlag > 0) IConst = Surface(SurfNum)%StormWinShadedConstruction
!END IF
!TotLay = Construct(IConst)%TotLayers
TotLay = Construct(ConstrNum)%TotLayers
IGap = 0
!****************************************************************************************************
!
! Inside and outside gas coefficients
!
!****************************************************************************************************
iprop(1, 1) = 1 ! air on outdoor side
frct(1, 1) = 1.0d0 ! pure air on outdoor side
nmix(1) = 1 ! pure air on outdoor side
iprop(nlayer + 1, 1) = 1 ! air on indoor side
frct(nlayer + 1, 1) = 1.0d0 ! pure air on indoor side
nmix(nlayer + 1) = 1 ! pure air on indoor side
!Simon: feed gas coefficients with air. This is necessary for tarcog because it is used on indoor and outdoor sides
GasType = GasCoeffsAir
wght(iprop(1, 1)) = GasWght(GasType)
gama(iprop(1, 1)) = GasSpecificHeatRatio(GasType)
DO ICoeff = 1,3
gcon(iprop(1, 1), ICoeff) = GasCoeffsCon(GasType,ICoeff)
gvis(iprop(1, 1), ICoeff) = GasCoeffsVis(GasType,ICoeff)
gcp(iprop(1, 1), ICoeff) = GasCoeffsCp (GasType,ICoeff)
END DO
! Fill window layer properties needed for window layer heat balance calculation
IGlass = 0
IGap = 0
DO Lay = 1,TotLay
LayPtr = Construct(ConstrNum)%LayerPoint(Lay)
IF(( Material(LayPtr)%Group == WindowGlass) .OR. (Material(LayPtr)%Group == WindowSimpleGlazing) ) THEN
IGlass = IGlass + 1
LayerType(IGlass) = 0 ! this marks specular layer type
thick(IGlass) = Material(LayPtr)%Thickness
scon(IGlass) = Material(LayPtr)%Conductivity
emis(2*IGlass-1) = Material(LayPtr)%AbsorpThermalFront
emis(2*IGlass) = Material(LayPtr)%AbsorpThermalBack
tir(2*IGlass-1) = Material(LayPtr)%TransThermal
tir(2*IGlass) = Material(LayPtr)%TransThermal
YoungsMod(IGlass) = Material(LayPtr)%YoungModulus
PoissonsRat(IGlass) = Material(LayPtr)%PoissonsRatio
else if (Material(LayPtr)%Group == ComplexWindowShade) then
if (CalcCondition == noCondition) then
if (Lay == 1) SurfaceWindow(SurfNum)%ShadingFlag = ExtShadeOn
if (Lay == TotLay) SurfaceWindow(SurfNum)%ShadingFlag = IntShadeOn
end if
IGlass = IGlass + 1
TempInt = Material(LayPtr)%ComplexShadePtr
LayerType(IGlass) = ComplexShade(TempInt)%LayerType
thick(IGlass) = ComplexShade(TempInt)%Thickness
scon(IGlass) = ComplexShade(TempInt)%Conductivity
emis(2*IGlass-1) = ComplexShade(TempInt)%FrontEmissivity
emis(2*IGlass) = ComplexShade(TempInt)%BackEmissivity
tir(2*IGlass-1) = ComplexShade(TempInt)%IRTransmittance
tir(2*IGlass) = ComplexShade(TempInt)%IRTransmittance
! This needs to be converted into correct areas. That can be done only after loading complete window data
Atop(IGlass) = ComplexShade(TempInt)%TopOpeningMultiplier
Abot(IGlass) = ComplexShade(TempInt)%BottomOpeningMultiplier
Al(IGlass) = ComplexShade(TempInt)%LeftOpeningMultiplier
Ar(IGlass) = ComplexShade(TempInt)%RightOpeningMultiplier
Ah(IGlass) = ComplexShade(TempInt)%FrontOpeningMultiplier
SlatThick(IGlass) = ComplexShade(TempInt)%SlatThickness
SlatWidth(IGlass) = ComplexShade(TempInt)%SlatWidth
SlatAngle(IGlass) = ComplexShade(TempInt)%SlatAngle
SlatCond(IGlass) = ComplexShade(TempInt)%SlatConductivity
SlatSpacing(IGlass) = ComplexShade(TempInt)%SlatSpacing
SlatCurve(IGlass) = ComplexShade(TempInt)%SlatCurve
else if(Material(LayPtr)%Group == ComplexWindowGap) THEN
IGap = IGap + 1
gap(IGap) = Material(LayPtr)%Thickness
presure(IGap) = Material(LayPtr)%Pressure
DeflectionPtr = Material(LayPtr)%DeflectionStatePtr
if(DeflectionPtr /= 0) then
GapDefMax(IGap) = DeflectionState(DeflectionPtr)%DeflectedThickness
else
GapDefMax(IGap) = gap(IGap)
end if
PillarPtr = Material(LayPtr)%SupportPillarPtr
if (PillarPtr.ne.0) then
SupportPlr(IGap) = 1
PillarSpacing(IGap) = SupportPillar(PillarPtr)%Spacing
PillarRadius(IGap) = SupportPillar(PillarPtr)%Radius
end if
GasPointer = Material(LayPtr)%GasPointer
nmix(IGap+1) = Material(GasPointer)%NumberOfGasesInMixture
DO IMix = 1,nmix(IGap+1)
!iprop(IGap+1, IMix) = Material(LayPtr)%GasType(IMix)
!iprop(IGap+1, IMix) = GetGasIndex(Material(LayPtr)%GasWGHT(IMix))
frct(IGap+1,IMix) = Material(GasPointer)%GasFract(IMix)
!Now has to build-up gas coefficients arrays. All used gasses should be stored into these arrays and
!to be correctly referenced by gap arrays
!First check if gas coefficients are already part of array. Duplicates are not necessary
CALL CheckGasCoefs(Material(GasPointer)%GasWGHT(IMix), iprop(IGap+1, IMix), wght, feedData)
IF (feedData) THEN
wght(iprop(IGap+1, IMix)) = Material(GasPointer)%GasWGHT(IMix)
gama(iprop(IGap+1, IMix)) = Material(GasPointer)%GasSpecHeatRatio(IMix)
DO i = 1, 3
gcon(iprop(IGap+1, IMix), i) = Material(GasPointer)%GasCON(IMix, i)
gvis(iprop(IGap+1, IMix), i) = Material(GasPointer)%GasVIS(IMix, i)
gcp(iprop(IGap+1, IMix), i) = Material(GasPointer)%GasCP(IMix, i)
END DO
END IF !IF feedData THEN
END DO
else
CALL ShowContinueError('Illegal layer type in Construction:ComplexFenestrationState.')
CALL ShowContinueError('Allowed object are:')
CALL ShowContinueError(' - WindowMaterial:Glazing')
CALL ShowContinueError(' - WindowMaterial:ComplexShade')
CALL ShowContinueError(' - WindowMaterial:Gap')
CALL ShowFatalError('halting because of error in layer definition for Construction:ComplexFenestrationState')
end if
END DO ! End of loop over glass, gap and blind/shade layers in a window construction
IF (CalcCondition == noCondition) THEN
! now calculate correct areas for multipliers
DO Lay = 1, nlayer
IF (LayerType(Lay) /= 0) THEN ! Layer is shading
! before changing multipliers, need to determine which one is dominant gap width
IF (Lay == 1) THEN ! Exterior shading device
dominantGapWidth = gap(Lay)
ELSE IF (Lay == nlayer) THEN ! Interior shading device
dominantGapWidth = gap(Lay - 1)
ELSE ! In-between shading device
dominantGapWidth = MIN(gap(Lay-1), gap(Lay))
END IF
Atop(Lay) = Atop(Lay) * dominantGapWidth * width
Abot(Lay) = Abot(Lay) * dominantGapWidth * width
Al(Lay) = Al(Lay) * dominantGapWidth * height
Ar(Lay) = Ar(Lay) * dominantGapWidth * height
Ah(Lay) = Ah(Lay) * width * height
END IF
END DO
END IF
!ThermalMod = 0
Debug_mode = 0
CalcSHGC = 0
Window_ID = ConstrNum
!vector of absorbed solar energy fractions for each layer.
asol = 0.0D0
! direct solar radiation
if (CalcCondition == noCondition) then
ShadeFlag = SurfaceWindow(SurfNum)%ShadingFlag
dir = QRadSWOutIncident(SurfNum) + QS(Surface(SurfNum)%Zone) ! TODO, check , !
! currently using Exterior beam plus diffuse solar incident on surface
! plus zone short wave. CHECK
!if (dir.ne.0.0d0) then
do IGlass = 1, nlayer
!IF (dir > 0.0D0 ) THEN
asol(IGLASS) = QRadSWwinAbs(SurfNum, IGLASS)
!ELSE
! asol(IGLASS) = 0.0D0
!ENDIF
end do
!end if
! Add contribution of IR from zone internal gains (lights, equipment and people). This is absorbed in zone-side layer and it
! is assumed that nothing is transmitted through
asol(nlayer) = asol(nlayer) + QRadThermInAbs(SurfNum)
presure = OutBaroPress
! Instead of doing temperature guess get solution from previous iteration. That should be much better than guess
do k = 1,2*nlayer
theta(k) = SurfaceWindow(SurfNum)%ThetaFace(k)
end do
end if
! Standard conditions run (winter and summer)
if (CalcCondition == winterCondition) then
tind = 294.15d0
tout = 255.15d0
hcout = 26.0d0
wso = 5.5d0
dir = 0.0d0
else if (CalcCondition == summerCondition) then
tind = 297.15d0
tout = 305.15d0
hcout = 15.0d0
wso = 2.75d0
dir = 783.0d0
CalcSHGC = 1
end if
! Common condition data
if (CalcCondition /= noCondition) then
trmin = tind
outir = 0.0d0
tsky = tout
wsi = 0.0d0
fclr = 1.0d0
ibc(1) = 0
ibc(2) = 0
presure = 101325.0d0
iwd = 0 ! Windward wind direction
isky = 0
esky = 1.0d0
height = 1.0d0
heightt = 1.0d0
width = 1.0d0
tilt = 90.0d0
! Just to make initial quess different from absolute zero
theta = 273.15d0
end if
! call TARCOG
CALL TARCOG90(nlayer, iwd, tout, tind, trmin, wso, wsi, dir, outir, isky, tsky, esky, fclr, VacuumPressure, &
VacuumMaxGapThickness, CalcDeflection, Pa, Pini, Tini, gap, GapDefMax, thick, scon, YoungsMod, PoissonsRat, &
tir, emis, totsol, tilt, asol, height, heightt, width, &
presure, iprop, frct, gcon, gvis, gcp, wght, gama, nmix, &
SupportPlr, PillarSpacing, PillarRadius, &
theta, LayerDef, q, qv, ufactor, sc, hflux, hcin, hcout, hrin, hrout, hin, hout, hcgap, hrgap, shgc, nperr, &
tarcogErrorMessage, shgct, tamb, troom, ibc, Atop, Abot, Al, Ar, Ah, SlatThick, SlatWidth, SlatAngle,&
SlatCond, SlatSpacing, SlatCurve, vvent,tvent, LayerType, nslice, LaminateA, LaminateB,&
sumsol, hg, hr, hs, he, hi,Ra,Nu, standard, ThermalMod, Debug_mode, Debug_dir, Debug_file, Window_ID, IGU_ID,&
ShadeEmisRatioOut, ShadeEmisRatioIn, ShadeHcRatioOut, ShadeHcRatioIn,&
HcUnshadedOut, HcUnshadedIn, Keff, ShadeGapKeffConv, SDScalar, CalcSHGC, NumOfIterations)
! process results from TARCOG
IF ((nperr > 0).and.(nperr < 1000)) THEN ! process error signal from tarcog
Call ShowSevereError('Window tarcog returned an error')
tarcogErrorMessage = 'message = "'//tarcogErrorMessage//'"'
Call ShowContinueErrorTimeStamp(tarcogErrorMessage)
if (CalcCondition == noCondition) then
CALL showContinueError('surface name = '//Surface(SurfNum)%Name)
end if
Call ShowContinueError('construction name = ' //Construct(ConstrNum)%Name)
CALL ShowFatalError('halting because of error in tarcog')
ELSE IF (CalcCondition == winterCondition) THEN
NominalU(ConstrNum) = ufactor
ELSE IF (CalcCondition == summerCondition) THEN
!tempInt = SurfaceWindow(SurfNum)%ComplexFen%CurrentState
!tempReal = SurfaceWindow(SurfNum)%ComplexFen%State(tempInt)%WinDiffTrans
!Sum1 = 0.0d0
!Sum2 = 0.0d0
!do j = 1 , ComplexWind(SurfNum)%Geom%Inc%NBasis !Incident ray loop
! Sum2 = Sum2 + ComplexWind(SurfNum)%Geom%Inc%Lamda (j)
! do m = 1 , ComplexWind(SurfNum)%Geom%Trn%NBasis !Outgoing ray loop
! Sum1 =Sum1 + ComplexWind(SurfNum)%Geom%Inc%Lamda(j) * ComplexWind(SurfNum)%Geom%Trn%Lamda(m) * &
! & Construct(ConstrNum)%BSDFInput%SolFrtTrans ( j , m )
! end do !Outgoing ray loop
!end do !Incident ray loop
!if (Sum2 > 0 ) THEN
! tempReal = Sum1/Sum2
!else
! tempReal = 0.
! call ShowWarningError ('BSDF--Inc basis has zero projected solid angle')
!endif
Construct(ConstrNum)%SummerSHGC = shgc
!Construct(SurfNum)%VisTransNorm = SurfaceWindow(SurfNum)%ComplexFen%State(tempInt)%WinDiffVisTrans
ELSE IF (CalcCondition == noCondition) THEN ! expect converged results...
! Window heat balance solution has converged.
SurfaceWindow(SurfNum)%WindowCalcIterationsRep = NumOfIterations
HConvIn(SurfNum) = hcin
! For interior shade, add convective gain from glass/shade gap air flow to zone convective gain;
! For all cases, get total window heat gain for reporting. See CalcWinFrameAndDividerTemps for
! contribution of frame and divider.
SurfInsideTemp = theta(2*nlayer) - KelvinConv
SurfOutsideTemp = theta(1) - KelvinConv
SurfOutsideEmiss = emis(1)
IncidentSolar = Surface(SurfNum)%Area * QRadSWOutIncident(SurfNum)
IF(ShadeFlag == IntShadeOn .OR. ShadeFlag == IntBlindOn) THEN
! Interior shade or blind
ConvHeatFlowNatural = -qv(nlayer) * height * width
SurfaceWindow(SurfNum)%ConvHeatFlowNatural = ConvHeatFlowNatural
WinGapConvHtFlowRep(SurfNum) = ConvHeatFlowNatural
WinGapConvHtFlowRepEnergy(SurfNum) = WinGapConvHtFlowRep(SurfNum) * TimeStepZone * SecInHour
! Window heat gain from glazing and shade/blind to zone. Consists of transmitted solar, convection
! from air exiting gap, convection from zone-side of shade/blind, net IR to zone from shade and net IR to
! zone from the glass adjacent to the shade/blind (zero if shade/blind IR transmittance is zero).
! Following assumes glazed area = window area (i.e., dividers ignored) in calculating
! IR to zone from glass when interior shade/blind is present.
ShadeArea = Surface(SurfNum)%Area + SurfaceWindow(SurfNum)%DividerArea
sconsh = scon(ngllayer+1) / thick(ngllayer+1)
nglfacep = nglface + 2
CondHeatGainShade = ShadeArea * sconsh * (theta(nglfacep-1) - theta(nglfacep))
EpsShIR1 = emis(nglface+1)
EpsShIR2 = emis(nglface+2)
TauShIR = tir(nglface+1)
RhoShIR1 = MAX(0.d0,1.d0-TauShIR-EpsShIR1)
RhoShIR2 = MAX(0.d0,1.d0-TauShIR-EpsShIR2)
RhoGlIR2 = 1.d0-emis(2*ngllayer)
ShGlReflFacIR = 1.d0-RhoGlIR2*RhoShIR1
NetIRHeatGainShade = ShadeArea * &
EpsShIR2*(sigma*theta(nglfacep)**4 - rmir) + &
EpsShIR1*(sigma*theta(nglfacep-1)**4 - rmir)*RhoGlIR2*TauShIR/ShGlReflFacIR
NetIRHeatGainGlass = ShadeArea * &
(emis(2*ngllayer)*TauShIR/ShGlReflFacIR) * (sigma*theta(2*ngllayer)**4 - rmir)
ConvHeatGainFrZoneSideOfShade = ShadeArea * hcin*(theta(nglfacep) - tind)
WinHeatGain(SurfNum) = WinTransSolar(SurfNum) + ConvHeatFlowNatural + ConvHeatGainFrZoneSideOfShade + &
NetIRHeatGainGlass + NetIRHeatGainShade
! store components for reporting
WinGainConvGlazShadGapToZoneRep(SurfNum) = ConvHeatFlowNatural ! result is in [W]
WinGainConvShadeToZoneRep(SurfNum) = ConvHeatGainFrZoneSideOfShade
WinGainIRGlazToZoneRep(SurfNum) = NetIRHeatGainGlass
WinGainIRShadeToZoneRep(SurfNum) = NetIRHeatGainShade
ELSE
! Interior shade or blind not present; innermost layer is glass
CondHeatGainGlass = Surface(SurfNum)%Area * scon(nlayer)/thick(nlayer) * (theta(2*nlayer-1)-theta(2*nlayer))
NetIRHeatGainGlass = Surface(SurfNum)%Area * emis(2*nlayer)*(sigma*theta(2*nlayer)**4 - rmir)
ConvHeatGainFrZoneSideOfGlass = Surface(SurfNum)%Area * hcin*(theta(2*nlayer) - tind)
WinHeatGain(SurfNum) = WinTransSolar(SurfNum) + ConvHeatGainFrZoneSideOfGlass + NetIRHeatGainGlass
! store components for reporting
WinGainConvGlazToZoneRep(SurfNum) = ConvHeatGainFrZoneSideOfGlass
WinGainIRGlazToZoneRep(SurfNum) = NetIRHeatGainGlass
! need to report convective heat flow from the gap in case of exterior shade
IF(ShadeFlag == ExtShadeOn) THEN
ConvHeatFlowNatural = -qv(2) * height * width ! qv(1) is exterior environment
WinGapConvHtFlowRep(SurfNum) = ConvHeatFlowNatural
WinGapConvHtFlowRepEnergy(SurfNum) = WinGapConvHtFlowRep(SurfNum) * TimeStepZone * SecInHour
END IF
END IF
! Add convective heat gain from airflow window
! Note: effect of fan heat on gap outlet temperature is neglected since fan power (based
! on pressure drop through the gap) is extremely small
!WinGapConvHtFlowRep(SurfNum) = 0.0d0
!WinGapConvHtFlowRepEnergy(SurfNum) = 0.0d0
TotAirflowGap = SurfaceWindow(SurfNum)%AirFlowThisTS * Surface(SurfNum)%Width
TAirflowGapOutlet = KelvinConv ! TODO Need to calculate this
TAirflowGapOutletC = TAirflowGapOutlet-KelvinConv
SurfaceWindow(SurfNum)%TAirflowGapOutlet = TAirflowGapOutletC
IF(SurfaceWindow(SurfNum)%AirFlowThisTS > 0.0d0) THEN
ConvHeatFlowForced = sum(qv) ! TODO. figure forced ventilation heat flow in Watts
WinGapConvHtFlowRep(SurfNum) = ConvHeatFlowForced
WinGapConvHtFlowRepEnergy(SurfNum) = WinGapConvHtFlowRep(SurfNum) * TimeStepZone * SecInHour
! Add heat from gap airflow to zone air if destination is inside air; save the heat gain to return
! air in case it needs to be sent to the zone (due to no return air determined in HVAC simulation)
IF(SurfaceWindow(SurfNum)%AirFlowDestination == AirFlowWindow_Destination_IndoorAir .or. &
SurfaceWindow(SurfNum)%AirFlowDestination == AirFlowWindow_Destination_ReturnAir) THEN
IF (SurfaceWindow(SurfNum)%AirflowSource == AirFlowWindow_Source_IndoorAir) THEN
InletAirHumRat = ZoneAirHumRat(ZoneNum)
ELSE ! AirflowSource = outside air
InletAirHumRat = OutHumRat
END IF
ZoneTemp = MAT(ZoneNum) ! this should be Tin (account for different reference temps)
CpAirOutlet = PsyCpAirFnWTdb(InletAirHumRat,TAirflowGapOutletC)
CpAirZone = PsyCpAirFnWTdb(ZoneAirHumRat(ZoneNum),ZoneTemp)
ConvHeatGainToZoneAir = TotAirflowGap * (CpAirOutlet*(TAirflowGapOutletC) - CpAirZone*ZoneTemp)
IF (SurfaceWindow(SurfNum)%AirFlowDestination == AirFlowWindow_Destination_IndoorAir) THEN
SurfaceWindow(SurfNum)%ConvHeatGainToZoneAir = ConvHeatGainToZoneAir
WinHeatGain(SurfNum) = WinHeatGain(SurfNum) + ConvHeatGainToZoneAir
ELSE
SurfaceWindow(SurfNum)%RetHeatGainToZoneAir = ConvHeatGainToZoneAir
END IF
END IF
! For AirflowDestination = ReturnAir in a controlled (i.e., conditioned) zone with return air, see CalcZoneLeavingConditions
! for calculation of modification of return-air temperature due to airflow from window gaps into return air.
END IF
! Correct WinHeatGain for interior diffuse shortwave (solar and shortwave from lights) transmitted
! back out window
ConstrNum = Surface(SurfNum)%Construction
!ConstrNumSh = Surface(SurfNum)%ShadedConstruction
!IF(SurfaceWindow(SurfNum)%StormWinFlag==1) THEN
! ConstrNum = Surface(SurfNum)%StormWinConstruction
! ConstrNumSh = Surface(SurfNum)%StormWinShadedConstruction
!END IF
!IF(ShadeFlag <= 0) THEN
TransDiff = Construct(ConstrNum)%TransDiff
!ELSE IF(ShadeFlag==IntShadeOn .OR. ShadeFlag==ExtShadeOn) THEN
! TransDiff = Construct(ConstrNum)%TransDiff
!ELSE IF(ShadeFlag==IntBlindOn .OR. ShadeFlag==ExtBlindOn .OR.ShadeFlag==BGBlindOn) THEN
! TransDiff = InterpSlatAng(SurfaceWindow(SurfNum)%SlatAngThisTS,SurfaceWindow(SurfNum)%MovableSlats, &
! Construct(ConstrNumSh)%BlTransDiff)
!ELSE IF(ShadeFlag == SwitchableGlazing) THEN
! TransDiff = InterpSW(SurfaceWindow(SurfNum)%SwitchingFactor,Construct(ConstrNum)%TransDiff, &
! Construct(ConstrNumSh)%TransDiff)
!END IF
WinHeatGain(SurfNum) = WinHeatGain(SurfNum) - QS(Surface(SurfNum)%Zone) * Surface(SurfNum)%Area * TransDiff
WinLossSWZoneToOutWinRep(SurfNum) = QS(Surface(SurfNum)%Zone) * Surface(SurfNum)%Area * TransDiff
IF(ShadeFlag==IntShadeOn.OR.ShadeFlag==ExtShadeOn) THEN
WinShadingAbsorbedSolar(SurfNum) = (SurfaceWindow(SurfNum)%ExtBeamAbsByShade + &
SurfaceWindow(SurfNum)%ExtDiffAbsByShade) * &
(Surface(SurfNum)%Area+SurfaceWindow(SurfNum)%DividerArea)
WinShadingAbsorbedSolarEnergy(SurfNum) = WinShadingAbsorbedSolar(SurfNum) * TimeStepZone * SecInHour
END IF
IF(SunIsUp) THEN
WinSysSolTransmittance(SurfNum) = WinTransSolar(SurfNum) / &
(QRadSWOutIncident(SurfNum)*(Surface(SurfNum)%Area+SurfaceWindow(SurfNum)%DividerArea)+0.0001d0)
WinSysSolAbsorptance(SurfNum) = (QRadSWwinAbsTot(SurfNum)+WinShadingAbsorbedSolar(SurfNum)) / &
(QRadSWOutIncident(SurfNum)*(Surface(SurfNum)%Area+SurfaceWindow(SurfNum)%DividerArea)+0.0001d0)
WinSysSolReflectance(SurfNum) = 1.0d0 - WinSysSolTransmittance(SurfNum) - WinSysSolAbsorptance(SurfNum)
ELSE
WinSysSolTransmittance(SurfNum) = 0.0d0
WinSysSolAbsorptance(SurfNum) = 0.0d0
WinSysSolReflectance(SurfNum) = 0.0d0
END IF
! Save hcv for use in divider calc with interior or exterior shade (see CalcWinFrameAndDividerTemps)
IF(ShadeFlag==IntShadeOn) SurfaceWindow(SurfNum)%ConvCoeffWithShade = 0.0d0
IF(ShadeFlag == IntShadeOn) THEN
SurfInsideTemp = theta(2*ngllayer+2) - KelvinConv
! Get properties of inside shading layer
ShadingLayPtr = Construct(ConstrNum)%LayerPoint(TotLay)
ShadingLayPtr = Material(ShadingLayPtr)%ComplexShadePtr
TauShadeIR = ComplexShade(ShadingLayPtr)%IRTransmittance
EpsShadeIR = ComplexShade(ShadingLayPtr)%FrontEmissivity
RhoShadeIR = MAX(0.d0,1.d0-TauShadeIR-EpsShadeIR)
! Get properties of glass next to inside shading layer
GlassLayPtr = Construct(ConstrNum)%LayerPoint(TotLay - 2)
EpsGlassIR = Material(GlassLayPtr)%AbsorpThermalBack
RhoGlassIR = 1 - EpsGlassIR
EffShBlEmiss = EpsShadeIR*(1.d0+RhoGlassIR*TauShadeIR/(1.d0-RhoGlassIR*RhoShadeIR))
SurfaceWindow(SurfNum)%EffShBlindEmiss = EffShBlEmiss
EffGlEmiss = EpsGlassIR*TauShadeIR/(1.d0-RhoGlassIR*RhoShadeIR)
SurfaceWindow(SurfNum)%EffGlassEmiss = EffGlEmiss
! EffShBlEmiss = InterpSlatAng(SurfaceWindow(SurfNum)%SlatAngThisTS,SurfaceWindow(SurfNum)%MovableSlats, &
! SurfaceWindow(SurfNum)%EffShBlindEmiss)
! EffGlEmiss = InterpSlatAng(SurfaceWindow(SurfNum)%SlatAngThisTS,SurfaceWindow(SurfNum)%MovableSlats, &
! SurfaceWindow(SurfNum)%EffGlassEmiss)
SurfaceWindow(SurfNum)%EffInsSurfTemp = (EffShBlEmiss * SurfInsideTemp + EffGlEmiss * (theta(2*ngllayer)-KelvinConv)) / &
(EffShBlEmiss + EffGlEmiss)
!ELSE
! SurfInsideTemp = theta(2*ngllayer) - TKelvin
!END IF
!IF(ShadeFlag == ExtShadeOn .OR. ShadeFlag == ExtBlindOn .OR. ShadeFlag == ExtScreenOn) THEN
!SurfOutsideTemp = theta(2*ngllayer+1) - TKelvin !this looks wrong.
ELSE
SurfOutsideTemp = theta(1) - KelvinConv
END IF
DO k = 1, nlayer
SurfaceWindow(SurfNum)%ThetaFace(2*k-1) = theta(2*k-1)
SurfaceWindow(SurfNum)%ThetaFace(2*k) = theta(2*k)
! temperatures for reporting
FenLaySurfTempFront(SurfNum, k) = theta(2*k-1) - KelvinConv
FenLaySurfTempBack(SurfNum, k) = theta(2*k) - KelvinConv
!thetas(k) = theta(k)
END DO
ENDIF
RETURN
END SUBROUTINE CalcComplexWindowThermal