Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer, | intent(in) | :: | nlayer | |||
integer, | intent(in) | :: | iwd | |||
real(kind=r64), | intent(in) | :: | tout | |||
real(kind=r64), | intent(in) | :: | tind | |||
real(kind=r64), | intent(in) | :: | wso | |||
real(kind=r64), | intent(in) | :: | wsi | |||
real(kind=r64), | intent(in) | :: | VacuumPressure | |||
real(kind=r64), | intent(in) | :: | VacuumMaxGapThickness | |||
real(kind=r64), | intent(inout) | :: | Ebsky | |||
real(kind=r64), | intent(out) | :: | tamb | |||
real(kind=r64), | intent(inout) | :: | Ebroom | |||
real(kind=r64), | intent(out) | :: | troom | |||
real(kind=r64), | intent(in), | dimension(MaxGap) | :: | gap | ||
real(kind=r64), | intent(in) | :: | height | |||
real(kind=r64), | intent(in) | :: | heightt | |||
real(kind=r64), | intent(in), | dimension(maxlay) | :: | scon | ||
real(kind=r64), | intent(in) | :: | tilt | |||
real(kind=r64), | intent(inout), | dimension(maxlay2) | :: | theta | ||
real(kind=r64), | intent(in), | dimension(maxlay1) | :: | Tgap | ||
real(kind=r64), | intent(inout), | dimension(maxlay2) | :: | Radiation | ||
real(kind=r64), | intent(in) | :: | Trmout | |||
real(kind=r64), | intent(in) | :: | Trmin | |||
integer, | intent(in), | dimension(maxlay1, maxgas) | :: | iprop | ||
real(kind=r64), | intent(in), | dimension(maxlay1, maxgas) | :: | frct | ||
real(kind=r64), | intent(in), | dimension(maxlay1) | :: | presure | ||
integer, | intent(in), | dimension(maxlay1) | :: | nmix | ||
real(kind=r64), | intent(in), | dimension(maxgas) | :: | wght | ||
real(kind=r64), | intent(in), | dimension(maxgas, 3) | :: | gcon | ||
real(kind=r64), | intent(in), | dimension(maxgas, 3) | :: | gvis | ||
real(kind=r64), | intent(in), | dimension(maxgas, 3) | :: | gcp | ||
real(kind=r64), | intent(in), | dimension(maxgas) | :: | gama | ||
integer, | intent(in) | :: | SupportPillar(maxlay) | |||
real(kind=r64), | intent(in), | dimension(maxlay) | :: | PillarSpacing | ||
real(kind=r64), | intent(in), | dimension(maxlay) | :: | PillarRadius | ||
real(kind=r64), | intent(out), | dimension(maxlay1) | :: | hgas | ||
real(kind=r64), | intent(out), | dimension(maxlay1) | :: | hcgas | ||
real(kind=r64), | intent(out), | dimension(maxlay1) | :: | hrgas | ||
real(kind=r64), | intent(inout) | :: | hcin | |||
real(kind=r64), | intent(inout) | :: | hcout | |||
real(kind=r64), | intent(in) | :: | hin | |||
real(kind=r64), | intent(in) | :: | hout | |||
integer, | intent(in) | :: | index | |||
integer, | intent(in), | dimension(2) | :: | ibc | ||
integer, | intent(inout) | :: | nperr | |||
character(len=*), | intent(inout) | :: | ErrorMessage | |||
real(kind=r64), | intent(inout) | :: | hrin | |||
real(kind=r64), | intent(inout) | :: | hrout | |||
real(kind=r64), | intent(inout), | dimension(maxlay) | :: | Ra | ||
real(kind=r64), | intent(inout), | dimension(maxlay) | :: | Nu |
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 hatter(nlayer, iwd, tout, tind, wso, wsi, VacuumPressure, VacuumMaxGapThickness, Ebsky, tamb, Ebroom, troom, &
gap, height, heightt, scon, tilt, theta, Tgap, Radiation, Trmout, Trmin, &
iprop, frct, presure, nmix, wght, gcon, gvis, gcp, gama, SupportPillar, PillarSpacing, PillarRadius, &
hgas, hcgas, hrgas, hcin, hcout, hin, hout, index, ibc, nperr, ErrorMessage, hrin, hrout, Ra, Nu )
!***********************************************************************
! This subroutine calculates the array of conductances/film coefficients used to model convection. The conductances/film
! coefficients are calculated as functions of temperature defined with the usual variable h and THEN are converted into an
! equivalent value interms of the black body emittance based on the surface
!***********************************************************************
! Inputs
! nlayer number of solid layers
! iwd wind direction
! tout outside temp in k
! tind inside temp in k
! wso wind speed in m/s
! wsi inside forced air speed m/s
! Ebsky ir flux from outside
! Ebroom ir flux from room
! Gout radiosity (ir flux) of the combined environment (sky+ground)
! Gin
! gap vector of gap widths in meters
! height IGU cavity height
! heightt
! thick glazing layer thickness
! scon Vector of conductivities of each glazing layer
! tilt Window tilt (in degrees)
! theta Vector of average temperatures
! Ebb
! Ebf
! iprop array of gap mixtures
! frct vector of mixture fractions
! presure
! hin Indoor Indoor combined film coefficient (if non-zero)
! hout Outdoor combined film coefficient (if non-zero)
! nmix vector of number of gasses in a mixture for each gap
! Ouputs
! hhat vector of all film coefficients (maxlay3)
! hgas vector of gap 'film' coeff.
! hcin Indoor convective surface heat transfer coefficient
! hcout Outdoor convective heat transfer coeff
! hrin Indoor radiative surface heat transfer coefficient
! hrout Outdoor radiative surface heat transfer coefficient
! hin Indoor combined film coefficient
! hout Outdoor combined film coefficient
! index iteration step
! ibc
!
! Inactives**
! wa - window azimuth (degrees, clockwise from south)
!
integer, intent(in) :: nlayer, index, iwd
integer, dimension(maxlay1, maxgas), intent(in) :: iprop
integer, dimension(maxlay1), intent(in) :: nmix
integer, dimension(2), intent(in) :: ibc
real(r64), dimension(maxgas), intent(in) :: wght, gama
real(r64), dimension(maxgas, 3), intent(in) :: gcon, gvis, gcp
real(r64), intent(in) :: tout, tind, wso, wsi, height, heightt, tilt, Trmout, Trmin, VacuumPressure, VacuumMaxGapThickness
real(r64), dimension(maxlay1, maxgas), intent(in) :: frct
real(r64), dimension(maxlay), intent(in) :: scon
real(r64), dimension(maxlay1), intent(in) :: presure
real(r64), dimension(maxlay2), intent(inout) :: theta, Radiation
real(r64), dimension(maxlay1), intent(in) :: Tgap
real(r64), intent(inout) :: Ebsky, Ebroom
real(r64), dimension(MaxGap), intent(in) :: gap
integer, intent(in) :: SupportPillar(maxlay)
real(r64), dimension(maxlay), intent(in) :: PillarSpacing, PillarRadius
real(r64), intent(in) :: hin, hout
real(r64), dimension(maxlay), intent(inout) :: Ra, Nu
real(r64), intent(inout) :: hrin, hrout, hcin, hcout
real(r64), dimension(maxlay1), intent(out):: hgas, hcgas, hrgas
real(r64), intent(out) :: troom, tamb
integer, intent(inout) :: nperr
character(len=*), intent(inout) :: ErrorMessage
integer :: i, k, nface
!character*(3) :: a
! common /props/ gcon(maxgas,3),gvis(maxgas,3),gcp(maxgas,3),grho(maxgas,3),wght(maxgas)
! evaluate convective/conductive components of gap grashof number, thermal conductivity and their derivatives:
nface = 2*nlayer
!do i = 1, nlayer
! j=2*i
! if ((Ebb(i)-Ebf(i)).eq.0) then
! theta(j) = theta(j) + tempCorrection
! Ebb(i) = sigma * (theta(j) ** 4)
! end if
! hhat(j) = scon(i)/thick(i) * (theta(j)-theta(j-1))/(Ebb(i)-Ebf(i))
!dr.....caluclate for laminate procedure
! if (nslice(i).gt.1) then
! if ((LaminateB(i).ne.0).and.((Ebb(i)-Ebf(i)).ne.0)) then
! hhat(j) = (theta(j)-theta(j-1))/(LaminateB(i) * (Ebb(i)-Ebf(i)))
! end if
! end if
!if (hhat(j).lt.0) then
! nperr = 6
! write(a, '(i3)') i
! ErrorMessage = 'Heat transfer coefficient based on emissive power in glazing layer is less than zero. Layer #'//trim(a)
! return
!end if
!end do
call filmg(tilt, theta, Tgap, nlayer, height, gap, iprop, frct, VacuumPressure, presure, nmix, &
wght, gcon, gvis, gcp, gama, hcgas, Ra, Nu, nperr, ErrorMessage)
if (.not.(GoAhead(nperr))) then
return
end if
!this is adding influence of pillar to hgas
call filmPillar(SupportPillar, scon, PillarSpacing, PillarRadius, nlayer, gap, hcgas, VacuumMaxGapThickness, &
nperr, ErrorMessage)
if (.not.(GoAhead(nperr))) then
return
end if
! adjust radiation coefficients
!hrgas = 0.0d0
do i=2, nlayer
k = 2*i-1
!if ((theta(k)-theta(k-1)) == 0) then
! theta(k-1) = theta(k-1) + tempCorrection
!end if
if ((theta(k)-theta(k-1)) /= 0) then
hrgas(i) = (Radiation(k) - Radiation(k-1))/(theta(k)-theta(k-1))
end if
hgas(i) = hcgas(i) + hrgas(i)
end do
! convective indoor film coeff:
if (ibc(2).le.0) then
call filmi(tind, theta(nface), nlayer, tilt, wsi, heightt, iprop, frct, presure, nmix, wght, gcon, gvis, gcp, hcin, ibc(2), &
nperr, ErrorMessage)
else if (ibc(2).eq.1) then
hcin = hin-hrin
!Simon: First iteration is with index = 0 and that means it should reenter iteration with whatever is provided as input
!else if (ibc(2).eq.2.and.index.eq.1) then
else if ((ibc(2).eq.2).and.(index.eq.0)) then
hcin=hin
end if
if (hcin.lt.0) then
nperr = 8
ErrorMessage = 'Hcin is less then zero.'
return
end if
hcgas(nlayer+1) = hcin
!hrin = 0.95*(Ebroom - Radiation(2*nlayer))/(Trmin-theta(2*nlayer))+0.05*hrin
hrin = (Ebroom - Radiation(2*nlayer))/(Trmin-theta(2*nlayer))
!if ((Theta(2*nlayer) - Trmin).ne.0) then
! hrin = sigma * emis(2*nlayer) * (Theta(2*nlayer)**4 - Trmin**4)/(Theta(2*nlayer) - Trmin)
!else
! Theta(2*nlayer) = Theta(2*nlayer) + tempCorrection
! hrin = sigma * emis(2*nlayer) * (Theta(2*nlayer)**4 - Trmin**4)/(Theta(2*nlayer) - Trmin)
!end if
hrgas(nlayer+1) = hrin
!hgas(nlayer+1) = hcgas(nlayer+1) + hrgas(nlayer+1)
troom = (hcin*tind+hrin*trmin)/(hcin+hrin)
! convective outdoor film coeff:
if (ibc(1).le.0) then
call film(tout,theta(1),wso,iwd,hcout,ibc(1))
else if (ibc(1).eq.1) then
hcout = hout-hrout
!Simon: First iteration is with index = 0 and that means it should reenter iteration with whatever is provided as input
!else if (ibc(1).eq.2.and.index.eq.1) then
else if ((ibc(1).eq.2).and.(index.eq.0)) then
hcout=hout
end if
if (hcout.lt.0) then
nperr = 9
ErrorMessage = 'Hcout is less than zero.'
return
end if
hcgas(1) = hcout
hrout = (Radiation(1) - Ebsky)/(theta(1)-Trmout)
!if ((Theta(1) - Trmout).ne.0) then
! hrout = sigma * emis(1) * (Theta(1)**4 - Trmout**4)/(Theta(1) - Trmout)
!else
! Theta(1) = Theta(1) + tempCorrection
! hrout = sigma * emis(1) * (Theta(1)**4 - Trmout**4)/(Theta(1) - Trmout)
!end if
hrgas(1) = hrout
!hgas(1) = hrout + hcout
tamb = (hcout*tout+hrout*trmout)/(hcout+hrout)
end subroutine hatter