Nodes of different colours represent the following:
Solid arrows point from a parent (sub)module to the submodule which is descended from it. Dashed arrows point from a module being used to the module or program unit using it. Where possible, edges connecting nodes are given different colours to make them easier to distinguish in large graphs.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer, | intent(in) | :: | Delt | |||
integer, | intent(in) | :: | I | |||
integer, | intent(in) | :: | Lay | |||
integer, | intent(in) | :: | Surf | |||
real(kind=r64), | intent(in), | DIMENSION(:) | :: | T | ||
real(kind=r64), | intent(inout), | DIMENSION(:) | :: | TT | ||
real(kind=r64), | intent(in), | DIMENSION(:) | :: | Rhov | ||
real(kind=r64), | intent(inout), | DIMENSION(:) | :: | RhoT | ||
real(kind=r64), | intent(inout), | DIMENSION(:) | :: | RH | ||
real(kind=r64), | intent(in), | DIMENSION(:) | :: | TD | ||
real(kind=r64), | intent(inout), | DIMENSION(:) | :: | TDT | ||
real(kind=r64), | intent(inout), | DIMENSION(:) | :: | EnthOld | ||
real(kind=r64), | intent(inout), | DIMENSION(:) | :: | EnthNew | ||
integer, | intent(in) | :: | TotNodes | |||
real(kind=r64), | intent(in) | :: | HMovInsul |
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 ExteriorBCEqns(Delt,I,Lay,Surf,T,TT,Rhov,RhoT,RH,TD,TDT,EnthOld,EnthNew,TotNodes,HMovInsul)
! SUBROUTINE INFORMATION:
! AUTHOR Richard Liesen
! DATE WRITTEN November, 2003
! MODIFIED B. Griffith 2010, fix adiabatic and other side surfaces
! May 2011, B. Griffith, P. Tabares
! November 2011 P. Tabares fixed problems with adiabatic walls/massless walls
! November 2011 P. Tabares fixed problems PCM stability problems
! RE-ENGINEERED Curtis Pedersen 2006
! PURPOSE OF THIS SUBROUTINE:
! <description>
! METHODOLOGY EMPLOYED:
! <description>
! REFERENCES:
! na
! USE STATEMENTS:
USE DataSurfaces, ONLY : OtherSideCondModeledExt, OSCM, HeatTransferModel_CondFD
USE DataHeatBalSurface , ONLY : QdotRadOutRepPerArea, QdotRadOutRep, QRadOutReport
IMPLICIT NONE ! Enforce explicit typing of all variables in this routine
! SUBROUTINE ARGUMENT DEFINITIONS:
INTEGER, INTENT(IN) :: Delt ! Time Increment
INTEGER, INTENT(IN) :: I ! Node Index
INTEGER, INTENT(IN) :: Lay ! Layer Number for Construction
INTEGER, INTENT(IN) :: Surf ! Surface number
REAL(r64),DIMENSION(:), INTENT(In) :: T !Old node Temperature in MFD finite difference solution
REAL(r64),DIMENSION(:), INTENT(InOut) :: TT !New node Temperature in MFD finite difference solution.
REAL(r64),DIMENSION(:), INTENT(In) :: Rhov !MFD Nodal Vapor Density[kg/m3] and is the old or last time step result.
REAL(r64),DIMENSION(:), INTENT(InOut) :: RhoT !MFD vapor density for the new time step.
REAL(r64),DIMENSION(:), INTENT(In) :: TD !The old dry Temperature at each node for the CondFD algorithm..
REAL(r64),DIMENSION(:), INTENT(InOut) :: TDT !The current or new Temperature at each node location for the CondFD solution..
REAL(r64),DIMENSION(:), INTENT(InOut) :: RH !Nodal relative humidity
REAL(r64),DIMENSION(:), INTENT(InOut) :: EnthOld ! Old Nodal enthalpy
REAL(r64),DIMENSION(:), INTENT(InOut) :: EnthNew ! New Nodal enthalpy
INTEGER, INTENT(IN) :: TotNodes ! Total nodes in layer
REAL(r64), INTENT(IN) :: HMovInsul ! Conductance of movable(transparent) insulation.
! SUBROUTINE PARAMETER DEFINITIONS:
! na
! INTERFACE BLOCK SPECIFICATIONS:
! na
! DERIVED TYPE DEFINITIONS:
! na
! SUBROUTINE LOCAL VARIABLE DECLARATIONS:
REAL(r64) :: QRadSWOutFD !Short wave radiation absorbed on outside of opaque surface
REAL(r64) :: Delx
INTEGER :: ConstrNum
INTEGER :: MatLay
INTEGER :: IndVarCol
INTEGER :: DepVarCol
REAL(r64) :: hconvo
REAL(r64) :: kt ! temperature dependent thermal conductivity, kt=ko +kt1(T-20)
REAL(r64) :: kto ! Base 20C thermal conductivity
REAL(r64) :: kt1 ! Thermal conductivity gradient coefficient where: kt=ko +kt1(T-20)
REAL(r64) :: Cpo ! Specific heat from idf
REAL(r64) :: Cp ! specific heat modified if PCM, otherwise equal to Cpo
REAL(r64) :: RhoS
REAL(r64) :: Toa
REAL(r64) :: Rhovo
REAL(r64) :: hmasso
REAL(r64) :: hgnd
REAL(r64) :: hrad
REAL(r64) :: hsky
REAL(r64) :: Tgnd
REAL(r64) :: Tsky
REAL(r64) :: Tia
REAL(r64) :: SigmaRLoc
REAL(r64) :: SigmaCLoc
REAL(r64) :: Rlayer
REAL(r64) :: QNetSurfFromOutside ! Combined outside surface net heat transfer terms
REAL(r64) :: TInsulOut ! Temperature of outisde face of Outside Insulation
REAL(r64) :: QRadSWOutMvInsulFD ! SW radiation at outside of Movable Insulation
INTEGER :: LayIn ! layer number for call to interior eqs
INTEGER :: NodeIn ! node number "I" for call to interior eqs
ConstrNum=Surface(surf)%Construction
!Boundary Conditions from Simulation for Exterior
hconvo = HConvExtFD(surf)
hmasso = HMassConvExtFD(surf)
hrad = HAirFD(surf)
hsky = HSkyFD(surf)
hgnd = HGrndFD(surf)
Toa = TempOutsideAirFD(surf)
rhovo = RhoVaporAirOut(surf)
Tgnd = TempOutsideAirFD(surf)
IF (Surface(Surf)%ExtBoundCond == OtherSideCondModeledExt) THEN
!CR8046 switch modeled rad temp for sky temp.
TSky = OSCM(Surface(Surf)%OSCMPtr)%TRad
QRadSWOutFD = 0.0D0 ! eliminate incident shortwave on underlying surface
ELSE
!Set the external conditions to local variables
QRadSWOutFD = QRadSWOutAbs(surf)
QRadSWOutMvInsulFD = QRadSWOutMvIns(Surf)
TSky = SkyTemp
ENDIF
Tia = Mat(Surface(Surf)%Zone)
SigmaRLoc=SigmaR(ConstrNum)
SigmaCLoc=SigmaC(ConstrNum)
MatLay = Construct(ConstrNum)%LayerPoint(Lay)
IF(Surface(Surf)%ExtBoundCond == Ground .or. IsRain)THEN
TDT(I) = Toa
TT(I) = Toa
RhoT(I)= rhovo
ELSEIF (Surface(Surf)%ExtBoundCond > 0) THEN
! this is actually the inside face of another surface, or maybe this same surface if adiabatic
! switch around arguments for the other surf and call routines as for interior side BC from opposite face
LayIn = Construct(Surface(Surface(Surf)%ExtBoundCond)%Construction)%TotLayers
nodeIn = ConstructFD(Surface(Surface(Surf)%ExtBoundCond)%Construction)%TotNodes + 1
IF (Surface(Surf)%ExtBoundCond == surf) THEN !adiabatic surface, PT addded since it is not the same as interzone wall
! as Outside Boundary Condition Object can be left blank.
CALL InteriorBCEqns(Delt,nodeIn,LayIn,Surf,SurfaceFD(Surf)%T, &
SurfaceFD(Surf)%TT, &
SurfaceFD(Surf)%Rhov, &
SurfaceFD(Surf)%RhoT, &
SurfaceFD(Surf)%RH, &
SurfaceFD(Surf)%TD, &
SurfaceFD(Surf)%TDT, &
SurfaceFD(Surf)%EnthOld, &
SurfaceFD(Surf)%EnthNew, &
SurfaceFD(Surf)%TDReport)
TDT(I) = SurfaceFD(Surf)%TDT(TotNodes + 1)
TT(I) = SurfaceFD(Surf)%TT(TotNodes + 1)
RhoT(I) = SurfaceFD(Surf)%RhoT(TotNodes + 1)
ELSE
CALL InteriorBCEqns(Delt,nodeIn,LayIn,Surface(Surf)%ExtBoundCond,SurfaceFD(Surface(Surf)%ExtBoundCond)%T, &
!potential-lkl-from old CALL InteriorBCEqns(Delt,nodeIn,LayIn,Surf,SurfaceFD(Surface(Surf)%ExtBoundCond)%T, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%TT, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%Rhov, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%RhoT, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%RH, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%TD, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%TDT, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%EnthOld, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%EnthNew, &
SurfaceFD(Surface(Surf)%ExtBoundCond)%TDReport)
TDT(I) = SurfaceFD(Surface(Surf)%ExtBoundCond)%TDT(TotNodes + 1)
TT(I) = SurfaceFD(Surface(Surf)%ExtBoundCond)%TT(TotNodes + 1)
RhoT(I) = SurfaceFD(Surface(Surf)%ExtBoundCond)%RhoT(TotNodes + 1)
ENDIF
! CALL InteriorBCEqns(Delt,nodeIn,Layin,Surface(Surf)%ExtBoundCond,SurfaceFD(Surf)%T, &
! SurfaceFD(Surf)%TT, &
! SurfaceFD(Surf)%Rhov, &
! SurfaceFD(Surf)%RhoT, &
! SurfaceFD(Surf)%RH, &
! SurfaceFD(Surf)%TD, &
! SurfaceFD(Surf)%TDT, &
! SurfaceFD(Surf)%EnthOld, &
! SurfaceFD(Surf)%EnthNew)
! now fill results from interior BC model eqns into local result for current call
! TDT(I) = SurfaceFD(Surface(Surf)%ExtBoundCond)%TDT(TotNodes + 1)
! TT(I) = SurfaceFD(Surface(Surf)%ExtBoundCond)%TT(TotNodes + 1)
! RhoT(I) = SurfaceFD(Surface(Surf)%ExtBoundCond)%RhoT(TotNodes + 1)
! TDT(I) = SurfaceFD(Surf)%TDT( i)
! TT(I) = SurfaceFD(Surf)%TT( i)
! RhoT(I) = SurfaceFD(Surf)%RhoT( i)
QNetSurfFromOutside = OpaqSurfInsFaceConductionFlux(Surface(Surf)%ExtBoundCond) !filled in InteriorBCEqns
! QFluxOutsideToOutSurf(Surf) = QnetSurfFromOutside
OpaqSurfOutsideFaceConductionFlux(Surf)= -QnetSurfFromOutside
OpaqSurfOutsideFaceConduction(Surf) = Surface(Surf)%Area * OpaqSurfOutsideFaceConductionFlux(Surf)
QHeatOutFlux(Surf) = QnetSurfFromOutside
ELSEIF (Surface(Surf)%ExtBoundCond <= 0) THEN ! regular outside conditions
!++++++++++++++++++++++++++++++++++++++++++++++++++++++
IF ( Surface(Surf)%HeatTransferAlgorithm == HeatTransferModel_CondFD ) THEN
! regular outside conditions
! Set Thermal Conductivity. Can be constant, simple linear temp dep or multiple linear segment temp function dep.
kto = Material(MatLay)%Conductivity ! 20C base conductivity
kt1 = MaterialFD(MatLay)%tk1 ! linear coefficient (normally zero)
kt= kto + kt1*((TDT(I)+TDT(I+1))/2.d0 - 20.d0)
IF( SUM(MaterialFD(MatLay)%TempCond(1:3,2)) >= 0.0d0) THEN ! Multiple Linear Segment Function
DepVarCol= 2 ! thermal conductivity
IndVarCol=1 !temperature
! Use average temp of surface and first node for k
kt =terpld(MaterialFD(MatLay)%numTempCond,MaterialFD(MatLay)%TempCond,(TDT(I)+TDT(I+1))/2.d0,IndVarCol,DepVarCol)
ENDIF
RhoS = Material(MatLay)%Density
Cpo = Material(MatLay)%SpecHeat
Cp=Cpo ! Will be changed if PCM
Delx = ConstructFD(ConstrNum)%Delx(Lay)
!Calculate the Dry Heat Conduction Equation
MatLay = Construct(ConstrNum)%LayerPoint(Lay)
IF (Material(MatLay)%ROnly .or. Material(MatLay)%Group == 1 ) THEN ! R Layer or Air Layer **********
! Use algebraic equation for TDT based on R
Rlayer = Material(MatLay)%Resistance
TDT(I)=(QRadSWOutFD*Rlayer + TDT(I+1) + hgnd*Rlayer*Tgnd + &
hconvo*Rlayer*Toa + hrad*Rlayer*Toa + &
hsky*Rlayer*Tsky)/ &
(1 + hconvo*Rlayer + hgnd*Rlayer + hrad*Rlayer + &
hsky*Rlayer)
IF ((TDT(I) > MaxSurfaceTempLimit) .OR. &
(TDT(I) < MinSurfaceTempLimit) ) THEN
TDT(I) = MAX(MinSurfaceTempLimit,MIN(MaxSurfaceTempLimit,TDT(I))) ! +++++ Limit Check
! CALL CheckFDSurfaceTempLimits(I,TDT(I))
ENDIF
ELSE ! Regular or phase change material layer
! check for phase change material
IF( SUM(MaterialFD(MatLay)%TempEnth(1:3,2)) >= 0.0d0) THEN ! phase change material, Use TempEnth Data to generate Cp
! CheckhT = Material(MatLay)%TempEnth ! debug
! Enthalpy function used to get average specific heat. Updated by GS so enthalpy function is followed.
DepVarCol= 2 ! enthalpy
IndVarCol=1 !temperature
EnthOld(I) =terpld(MaterialFD(MatLay)%numTempEnth,MaterialFD(MatLay)%TempEnth,TD(I),IndVarCol,DepVarCol)
EnthNew(I) =terpld(MaterialFD(MatLay)%numTempEnth,MaterialFD(MatLay)%TempEnth,TDT(I),IndVarCol,DepVarCol)
IF (EnthNew(I)==EnthOld(I)) THEN
Cp=Cpo
ELSE
Cp=MAX(Cpo,(EnthNew(I) -EnthOld(I))/(TDT(I)-TD(I)))
END IF
ELSE
Cp = Cpo
END IF ! Phase Change Material option
! Choose Regular or Transparent Insulation Case
IF(HMovInsul <= 0.0d0 ) THEN ! regular case
SELECT CASE (CondFDSchemeType)
CASE (CrankNicholsonSecondOrder)
! Second Order equation
TDT(I)= (1.0d0*QRadSWOutFD + (0.5d0*Cp*Delx*RhoS*TD(I))/DelT + (0.5d0*kt*(-1.0d0*TD(I) + TD(I+1)))/Delx &
+ (0.5d0*kt*TDT(I+1))/Delx + 0.5d0*hgnd*Tgnd + 0.5d0*hgnd*(-1.0d0*TD(I) + Tgnd) + 0.5d0*hconvo*Toa + &
0.5d0*hrad*Toa &
+ 0.5d0*hconvo*(-1.d0*TD(I) + Toa) + 0.5d0*hrad*(-1.d0*TD(I) + Toa) + 0.5d0*hsky*Tsky + &
0.5d0*hsky*(-1.d0*TD(I) + Tsky))/ &
(0.5d0*hconvo + 0.5d0*hgnd + 0.5d0*hrad + 0.5d0*hsky + (0.5d0*kt)/Delx + (0.5d0*Cp*Delx*RhoS)/DelT)
!feb2012 TDT(I)= (1.0d0*QRadSWOutFD + (0.5d0*Cp*Delx*RhoS*TD(I))/DelT + (0.5d0*kt*(-1.0d0*TD(I) + TD(I+1)))/Delx &
!feb2012 + (0.5d0*kt*TDT(I+1))/Delx + 0.5d0*hgnd*Tgnd + 0.5d0*hgnd*(-1.0d0*TD(I) + Tgnd) + 0.5d0*hconvo*Toa + &
!feb2012 0.5d0*hrad*Toa &
!feb2012 + 0.5d0*hconvo*(-1.d0*TD(I) + Toa) + 0.5d0*hrad*(-1.d0*TD(I) + Toa) + 0.5d0*hsky*Tsky + &
!feb2012 0.5d0*hsky*(-1.d0*TD(I) + Tsky))/ &
!feb2012 (0.5d0*hconvo + 0.5d0*hgnd + 0.5d0*hrad + 0.5d0*hsky + (0.5d0*kt)/Delx + (0.5d0*Cp*Delx*RhoS)/DelT)
CASE (FullyImplicitFirstOrder)
! First Order
TDT(I)=(2.0d0*Delt*Delx*QRadSWOutFD + Cp*Delx**2*RhoS*TD(I) &
+ 2.0d0*Delt*kt*TDT(I+1) + 2.0d0*Delt*Delx*hgnd*Tgnd &
+ 2.0d0*Delt*Delx*hconvo*Toa + 2.0d0*Delt*Delx*hrad*Toa + 2.0d0*Delt*Delx*hsky*Tsky)/ &
(2.0d0*Delt*Delx*hconvo + 2.0d0*Delt*Delx*hgnd + 2.0d0*Delt*Delx*hrad + &
2.0d0*Delt*Delx*hsky + 2.0d0*Delt*kt + Cp*Delx**2*RhoS)
END SELECT
ELSE IF ( HMovInsul > 0.0d0) Then ! Transparent insulation on outside
! Transparent insulaton additions
!Movable Insulation Layer Outside surface temp
TInsulOut = (QRadSWOutMvInsulFD + hgnd*Tgnd + HmovInsul*TDT(I) + &
hconvo*Toa + hrad*Toa + hsky*Tsky)/ &
(hconvo + hgnd + HmovInsul + hrad + hsky)
!List(List(Rule(TDT,(2*Delt*Delx*QradSWOutAbs +
!- Cp*Delx**2*Rhos*TD + 2*Delt*kt*TDTP1 +
!- 2*Delt*Delx*HmovInsul*Tiso)/
!- (2*Delt*Delx*HmovInsul + 2*Delt*kt + Cp*Delx**2*Rhos))))
! Wall first node temperature behind Movable insulation
SELECT CASE (CondFDSchemeType)
CASE (CrankNicholsonSecondOrder)
TDT(I)= (2*Delt*Delx*QradSWOutFD + &
Cp*Delx**2*Rhos*TD(I) + 2*Delt*kt*TDT(I+1) + &
2*Delt*Delx*HmovInsul*TInsulOut)/ &
(2*Delt*Delx*HmovInsul + 2*Delt*kt + Cp*Delx**2*Rhos)
CASE (FullyImplicitFirstOrder)
! Currently same as Crank Nicholson, need fully implicit formulation
TDT(I)= (2*Delt*Delx*QradSWOutFD + &
Cp*Delx**2*Rhos*TD(I) + 2*Delt*kt*TDT(I+1) + &
2*Delt*Delx*HmovInsul*TInsulOut)/ &
(2*Delt*Delx*HmovInsul + 2*Delt*kt + Cp*Delx**2*Rhos)
END SELECT
End IF ! Regular layer or Movable insulation cases
IF ((TDT(I) > MaxSurfaceTempLimit) .OR. &
(TDT(I) < MinSurfaceTempLimit) ) THEN
TDT(I) = Max(MinSurfaceTempLimit,Min(MaxSurfaceTempLimit,TDT(I))) ! +++++ Limit Check
! CALL CheckFDSurfaceTempLimits(I,TDT(I))
ENDIF
END IF ! R layer or Regular layer
END IF !End IF--ELSE SECTION (regular detailed FD part or SigmaR SigmaC part
! Determine net heat flux to ooutside face.
!One formulation that works for Fully Implicit and CrankNicholson and massless wall
QNetSurfFromOutside = QRadSWOutFD + (hgnd*(-TDT(I) + Tgnd) + &
hconvo*(-TDT(I) + Toa) + hrad*(-TDT(I) + Toa) + hsky*(-TDT(I) + Tsky))
!Same sign convention as CTFs
OpaqSurfOutsideFaceConductionFlux(Surf) = -QnetSurfFromOutside
OpaqSurfOutsideFaceConduction(Surf) = Surface(Surf)%Area * OpaqSurfOutsideFaceConductionFlux(Surf)
!Report all outside BC heat fluxes
QdotRadOutRepPerArea(Surf) = -(hgnd*(TDT(I) - Tgnd) + hrad*(TDT(I) - Toa) + hsky*( TDT(I) - Tsky))
QdotRadOutRep(Surf) = Surface(Surf)%Area * QdotRadOutRepPerArea(Surf)
QRadOutReport(Surf) = QdotRadOutRep(Surf) * SecInHour * TimeStepZone
END IF !End IF --ELSE SECTION (regular BC part of the ground and Rain check)
RETURN
END SUBROUTINE ExteriorBCEqns