Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer, | intent(in) | :: | SDLayerIndex | |||
integer, | intent(in), | dimension(2) | :: | ibc | ||
real(kind=r64), | intent(in) | :: | tout | |||
real(kind=r64), | intent(in) | :: | tind | |||
integer, | intent(in) | :: | nlayer | |||
real(kind=r64), | intent(in), | dimension(maxlay2) | :: | theta | ||
real(kind=r64), | intent(in) | :: | wso | |||
real(kind=r64), | intent(in) | :: | wsi | |||
integer, | intent(in) | :: | iwd | |||
real(kind=r64), | intent(in) | :: | height | |||
real(kind=r64), | intent(in) | :: | heightt | |||
real(kind=r64), | intent(in) | :: | tilt | |||
real(kind=r64), | intent(in), | dimension(maxlay) | :: | thick | ||
real(kind=r64), | intent(in), | dimension(MaxGap) | :: | gap | ||
real(kind=r64), | intent(in) | :: | hout | |||
real(kind=r64), | intent(in) | :: | hrout | |||
real(kind=r64), | intent(in) | :: | hin | |||
real(kind=r64), | intent(in) | :: | hrin | |||
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 | ||
integer, | intent(in) | :: | index | |||
real(kind=r64), | intent(in) | :: | SDScalar | |||
real(kind=r64), | intent(in), | dimension(maxlay) | :: | Ebf | ||
real(kind=r64), | intent(in), | dimension(maxlay) | :: | Ebb | ||
real(kind=r64), | intent(inout), | dimension(maxlay) | :: | hgas | ||
real(kind=r64), | intent(inout), | dimension(maxlay3) | :: | hhat | ||
integer, | intent(inout) | :: | nperr | |||
character(len=*), | intent(inout) | :: | ErrorMessage |
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 adjusthhat(SDLayerIndex, ibc, tout, tind, nlayer, theta, wso, wsi, iwd, height, heightt, tilt, &
thick, gap, hout, hrout, hin, hrin, iprop, frct, presure, nmix, wght, gcon, gvis, gcp, &
index, SDScalar, Ebf, Ebb, hgas, hhat, nperr, ErrorMessage)
!********************************************************************
! Modifies hhat, hgas coefficients around SD layers
!********************************************************************
integer, intent(in) :: SDLayerIndex, nlayer, iwd, index
integer, dimension(2), intent(in) :: ibc
integer, dimension(maxlay1, maxgas), intent(in) :: iprop
integer, dimension(maxlay1), intent(in) :: nmix
real(r64), dimension(maxgas), intent(in) :: wght
real(r64), dimension(maxgas, 3), intent(in) :: gcon, gvis, gcp
real(r64), intent(in) :: tout, tind, wso, wsi, height, heightt, tilt
real(r64), dimension(maxlay2), intent(in) :: theta
real(r64), dimension(maxlay), intent(in) :: thick
real(r64), dimension(MaxGap), intent(in) :: gap
real(r64), intent(in) :: hout, hrout, hin, hrin, SDScalar
real(r64), dimension(maxlay1, maxgas), intent(in) :: frct
real(r64), dimension(maxlay1), intent(in) :: presure
real(r64), dimension(maxlay), intent(in) :: Ebf, Ebb
real(r64), dimension(maxlay), intent(inout) :: hgas
real(r64), dimension(maxlay3), intent(inout) :: hhat
integer, intent(inout):: nperr
character(len=*), intent(inout) :: ErrorMessage
real(r64) :: hc_NOSD, hc_0, hc_1, hc_alpha, hhat_alpha
real(r64) :: hc_1_1, hc_1_2, hc_alpha1, hc_alpha2, hhat_alpha1, hhat_alpha2
real(r64), dimension(maxgas) :: frctg
integer :: i, j, k, l
integer, dimension(maxgas) :: ipropg
real(r64) :: tmean, con, visc, dens, cp, pr, delt, gap_NOSD, rayl, asp, gnu
!bi... Step 1: Calculate hc as if there was no SD here
if (SDLayerIndex.eq.1) then
!car SD is the first layer (outdoor)
! calc hc_0 as hcout:
! convective outdoor film coeff:
if (ibc(1).le.0) then
call film(tout,theta(3),wso,iwd,hc_NOSD,ibc(1))
else if (ibc(1).eq.1) then
hc_NOSD = hout-hrout
else if ((ibc(1).eq.2).and.(index.eq.1)) then
hc_NOSD=hout
end if
if (hc_NOSD.lt.0) then
nperr = 9
ErrorMessage = 'Hcout is out of range.'
return
end if
hc_0 = hc_NOSD*(theta(3)-tout)/(theta(3)-theta(2))
hc_1 = hgas(2)
hc_alpha = SDScalar*(hc_1 - hc_0) + hc_0
hhat_alpha = hc_alpha*(theta(3)-theta(2)) / (Ebf(2)-Ebb(1))
hgas(2) = hc_alpha
hhat(3) = hhat_alpha
else if (SDLayerIndex.eq.nlayer) then
!car SD is the last layer (indoor)
! calc hc_0 as hcin:
! convective indoor film coeff:
if (ibc(2).le.0) then
call filmi(tind, theta(2*nlayer-2), nlayer, tilt, wsi, heightt, iprop, frct, presure, nmix, wght, gcon, gvis, gcp, &
hc_NOSD, ibc(2), nperr, ErrorMessage)
else if (ibc(2).eq.1) then
hc_NOSD = hin-hrin
else if (ibc(2).eq.2.and.index.eq.1) then
hc_NOSD = hin
end if
if (hc_NOSD.lt.0) then
nperr = 8
ErrorMessage = 'Hcin is out of range.'
return
end if
! hgas(2*nlayer+1) = hcin
hc_0 = hc_NOSD*(tind-theta(2*nlayer-2))/(theta(2*nlayer-1)-theta(2*nlayer-2))
hc_1 = hgas(nlayer+1)
hc_alpha = SDScalar*(hc_1 - hc_0) + hc_0
hhat_alpha = hc_alpha*(theta(2*nlayer-1) - theta(2*nlayer-2)) / (Ebf(nlayer)-Ebb(nlayer-1))
hgas(nlayer+1) = hc_alpha
hhat(2*nlayer-1) = hhat_alpha
else
!car SD is in between glazing
! calc hc_NOSD as hcgas:
j = 2*SDLayerIndex-2
k = j+3
! determine the gas properties for this "gap":
tmean = (theta(j)+theta(k))/2.0d0
delt = abs(theta(j)-theta(k))
i = SDLayerIndex
do l=1, nmix(i+1)
ipropg(l) = iprop(i+1,l)
frctg(l) = frct(i+1,l)
end do
call gasses90(tmean, ipropg, frctg, presure(i+1), nmix(i+1), wght, gcon, gvis, gcp, con, visc, dens, cp, pr, &
ISO15099, nperr, ErrorMessage)
gap_NOSD = gap(SDLayerIndex-1) + gap(SDLayerIndex) + thick(SDLayerIndex)
! determine the Rayleigh number:
rayl = GravityConstant * gap_NOSD**3 * delt* cp * dens**2 / (tmean*visc*con)
asp = height / gap_NOSD
! determine the Nusselt number:
call nusselt(tilt, rayl, asp, gnu, nperr, ErrorMessage)
! calculate effective conductance of the gap
hc_NOSD = con/gap_NOSD*gnu
i = SDLayerIndex
j = 2 * i
!car changed hc interpolation with temperatures by inverting ratios
!car hc_0_1 = hc_NOSD*abs((theta(2*i-1) - theta(2*i-2))/(theta(2*i+1) - theta(2*i-2)))
!car hc_0_2 = hc_NOSD*abs((theta(2*i+1) - theta(2*i))/(theta(2*i+1) - theta(2*i-2)))
hc_1_1 = hgas(i)
hc_1_2 = hgas(i+1)
!car changed how SDscalar interpolation is done. Instead of hc_0_1 and hc_0_2, hc_NOSD is used
hc_alpha1 = SDScalar*(hc_1_1 - hc_NOSD) + hc_NOSD
hc_alpha2 = SDScalar*(hc_1_2 - hc_NOSD) + hc_NOSD
hhat_alpha1 = hc_alpha1*(theta(j-1) - theta(j-2)) / (Ebf(i)-Ebb(i-1))
hhat_alpha2 = hc_alpha2*(theta(j+1) - theta(j)) / (Ebf(i+1)-Ebb(i))
hgas(i) = hc_alpha1
hgas(i+1) = hc_alpha2
hhat(j-1) = hhat_alpha1
hhat(j+1) = hhat_alpha2
end if
end subroutine adjusthhat