Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real(kind=r64), | intent(in) | :: | tair | |||
real(kind=r64), | intent(in) | :: | t | |||
integer, | intent(in) | :: | nlayer | |||
real(kind=r64), | intent(in) | :: | tilt | |||
real(kind=r64), | intent(in) | :: | wsi | |||
real(kind=r64), | intent(in) | :: | height | |||
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(out) | :: | hcin | |||
integer, | intent(in) | :: | ibc | |||
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.
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 filmi(tair,t,nlayer,tilt,wsi,height,iprop,frct,presure,nmix,wght,gcon,gvis,gcp,hcin,ibc,nperr,ErrorMessage)
!***********************************************************************
! purpose to evaluate heat flux at indoor surface of window using still air correlations (Curcija and Goss 1993)
! found in SPC142 equations 5.43 - 5.48.
!***********************************************************************
! Input
! tair - room air temperature
! t - inside surface temperature
! nlayer number of glazing layers
! tilt - the tilt of the glazing in degrees
! wsi - room wind speed (m/s)
! height - window height
! iprop
! frct
! presure
! nmix vector of number of gasses in a mixture for each gap
! Output
! hcin - indoor convecive heat transfer coeff
! If there is forced air in the room than use SPC142 corelation 5.49 to calculate the room side film coefficient.
real(r64), intent(in) :: tair, t, tilt, wsi, height
real(r64), dimension(maxlay1, maxgas), intent(in) :: frct
real(r64), dimension(maxlay1), intent(in) :: presure
real(r64), dimension(maxgas), intent(in) :: wght
real(r64), dimension(maxgas, 3), intent(in) :: gcon, gvis, gcp
integer, intent(in) :: nlayer, ibc
integer, dimension(maxlay1, maxgas), intent(in) :: iprop
integer, dimension(maxlay1), intent(in) :: nmix
real(r64), intent(out) :: hcin
integer, intent(inout) :: nperr
character(len=*), intent(inout) :: ErrorMessage
real(r64), dimension(maxgas) :: frcti
integer :: j
integer, dimension(maxgas) :: ipropi
real(r64) :: tiltr, tmean, delt, con, visc, dens, cp, pr, gr, RaCrit, RaL, Gnui
if (wsi .gt. 0.0d0) then ! main IF
select case (ibc)
case (0)
hcin = 4.0d0 + 4.0d0 * wsi
case (-1)
hcin = 5.6d0 + 3.8d0 * wsi ! SPC142 correlation
return
end select
else ! main IF - else
tiltr = tilt*2.0d0*pi/360.0d0 ! convert tilt in degrees to radians
tmean = tair + 0.25d0 * (t - tair)
delt = ABS(tair-t)
do j=1, nmix(nlayer+1)
ipropi(j) = iprop(nlayer+1,j)
frcti(j) = frct(nlayer+1,j)
end do
call gasses90(tmean, ipropi, frcti, presure(nlayer+1), nmix(nlayer+1), wght, gcon, gvis, gcp, con, visc, dens, cp, pr, &
ISO15099, nperr, ErrorMessage)
! Calculate grashoff number:
! The grashoff number is the Rayleigh Number (equation 5.29) in SPC142 divided by the Prandtl Number (prand):
gr = GravityConstant * height**3 * delt*dens**2 / (tmean*visc**2)
RaL = gr*pr
! write(*,*)' RaCrit,RaL,gr,pr '
! write(*,*) RaCrit,RaL,gr,pr
if ((0.0d0.le.tilt).and.(tilt.lt.15.0d0)) then ! IF no. 1
Gnui = 0.13d0 * RaL**(1.0d0/3.0d0)
else if ((15.0d0.le.tilt).and.(tilt.le.90.0d0)) then
! if the room air is still THEN use equations 5.43 - 5.48:
RaCrit = 2.5d5 * (EXP(0.72d0 * tilt) / SIN(tiltr))**0.2d0
if (RaL.le.RaCrit) THEN ! IF no. 2
Gnui = 0.56d0 * (RaL * SIN(tiltr))**0.25d0
! write(*,*) ' Nu ', Gnui
else
!Gnui = 0.13*(RaL**0.3333 - RaCrit**0.3333) + 0.56*(RaCrit*sin(tiltr))**0.25
Gnui = 0.13d0 * (RaL**(1.0d0/3.0d0) - RaCrit**(1.0d0/3.0d0)) + 0.56d0 * (RaCrit*sin(tiltr))**0.25d0
end if ! end if no. 2
else if ((90.0d0.lt.tilt).and.(tilt.le.179.0d0)) then
Gnui = 0.56d0 * (RaL*SIN(tiltr))**0.25d0
else if ((179.0d0.lt.tilt).and.(tilt.le.180.0d0)) then
Gnui = 0.58d0 * RaL**(1/3.0d0)
end if ! end if no. 1
! write(*,*) ' RaL ', RaL, ' RaCrit', RaCrit
! write(*,*)' Nusselt Number ',Gnui
hcin = Gnui * (con/height)
! hin = 1.77*(abs(t-tair))**0.25
end if ! end main IF
end subroutine filmi