Inside Heat Balance[LINK]
The heart of the heat balance method is the internal heat balance involving the inside faces of the zone surfaces. This heat balance is generally modeled with four coupled heat transfer components: 1) conduction through the building element, 2) convection to the air, 3) short wave radiation absorption and reflectance and 4) longwave radiant interchange. The incident short wave radiation is from the solar radiation entering the zone through windows and emittance from internal sources such as lights. The longwave radiation interchange includes the absorption and emittance of low temperature radiation sources, such as all other zone surfaces, equipment, and people.
The heat balance on the inside face can be written as follows:
q′′LWX+q′′SW+q′′LWS+q′′ki+q′′sol+q′′conv=0
where:
q′′LWX = Net longwave radiant exchange flux between zone surfaces.
q′′SW = Net short wave radiation flux to surface from lights.
q′′LWS = Longwave radiation flux from equipment in zone.
q′′ki = Conduction flux through the wall.
q′′sol = Transmitted solar radiation flux absorbed at surface.
q′′conv = Convective heat flux to zone air.
Each of these heat balance components is introduced briefly below.
Internal LongWave Radiation Exchange[LINK]
LW Radiation Exchange Among Zone Surfaces[LINK]
There are two limiting cases for internal LW radiation exchange that are easily modeled:
The limiting case of completely absorbing air has been used for load calculations and also in some energy analysis calculations. This model is attractive because it can be formulated simply using a combined radiation and convection heat transfer coefficient from each surface to the zone air. However, it oversimplifies the zone surface exchange problem, and as a result, the heat balance formulation in EnergyPlus treats air as completely transparent. This means that it does not participate in the LW radiation exchange among the surfaces in the zone. The model, which considers room air to be completely transparent, is reasonable physically because of the low water vapor concentrations and the short mean path lengths. It also permits separating the radiant and convective parts of the heat transfer at the surface, which is an important attribute of the heat balance method.
EnergyPlus uses a grey interchange model for the longwave radiation among zone surfaces. This model is based on the “ScriptF” concept developed by Hottel (Hottel and Sarofim, Radiative Transfer, Chapter 3, McGraw Hill, 1967). This procedure relies on a matrix of exchange coefficients between pairs of surfaces that include all exchange paths between the surfaces. In other words all reflections, absorptions and reemissions from other surfaces in the enclosure are included in the exchange coefficient, which is called ScriptF. The major assumptions are that all surface radiation properties are grey and all radiation is diffuse. Both assumptions are reasonable for building zone interchange.
The ScriptF coefficients are developed by starting with the traditional direct radiation view factors. In the case of building rooms and zones, there are several complicating factors in finding the direct view factors—the main one being that the location of surfaces such as thermal mass representing furniture and partitions are not known. The other limitation is that the exact calculation of direct view factors is computationally very intensive even if the positions of all surfaces are known. Accordingly, EnergyPlus uses a procedure to approximate the direct view factors. The procedure has two steps:
1) Determine the total area of other surfaces “seen” by a surface.
2) Approximate the direct view factor from surface 1 to surface 2 as the ratio of the area of surface 2 to the total area “seen” by surface 1.
The determination of the “seen” area has several constraints:
No surface sees itself.
All surfaces see thermal mass surfaces.
No surface facing within 10 degrees of another surface is seen by the other surface.
All surfaces see roofs, floors and ceilings (subject to the preceding facing direction constraint).
Because the approximate view factors may not satisfy the basic requirements of reciprocity (two surfaces should exchange equal amounts of heat in each direction), and completeness (every surface should have a direct view factor sum of 1.0), EnergyPlus does a view factor fix operation before they are used in the ScriptF determination. Normally both of the requirements are satisfied, but in some special situations they are not, and special rules are applied.
If a user includes less than four surfaces in a zone, the enforcement of “completeness” is not strictly possible because it is not possible to have a completed enclosure with less than four surfaces. In this situation, reciprocity is enforced initially and then view factors from each surface are summed. When dealing with this few surfaces, the result of using approximate view factors is that once reciprocity is enforced that the largest surface in the zone will end up having view factors that sum up to greater than unity. This is not physically possible. So, to correct this and maintain completeness for one surface, the entire view factor matrix is divided by the largest summation of view factors for any surface in the zone. By dividing the entire view factor matrix by the same value, the reciprocity established by the standard algorithm is maintained and completeness for one surface is established. All of the other surfaces will not have view factors that sum up to unity and thus will not satisfy completeness. However, these other surfaces will not have view factors that sum up to greater than unity thus avoiding a physically impossible situation.
If the area of one surface in a zone is greater than the sum of the areas of all other surfaces, reciprocity only is enforced, but sometimes, for very large surfaces, that enforcement becomes impossible, and the view factors are modified so that only the large surface is seen by very small surfaces.
Warning messages are produced for both of these cases, and the results should be examined very carefully to ascertain that they are reasonable. The suggested action for the second case (the extralarge surface) is to divide the large surface into several smaller surfaces; then the enclosure will be treated as normal.
Once the ScriptF coefficients are determined, the longwave radiant exchange is calculated for each surface using:
qi,j=AiFi,j(T4i−T4j)
where Fi,j is the ScriptF between surfaces i and j.
Thermal Mass and Furniture[LINK]
Furniture in a zone has the effect of increasing the amount of surface area that can participate in the radiation and convection heat exchanges. It also adds participating thermal mass to the zone. These two changes both affect the response to temperature changes in the zone and also affect the heat extraction characteristics.
The proper modeling of furniture is an area that needs further research, but the heat balance formulation allows the effect to be modeled in a realistic manner by including the furniture surface area and thermal mass in the heat exchange process.
LW Radiation From Internal Sources[LINK]
The traditional model for this source is to define a radiative/convective split for the heat introduced into a zone from equipment. The radiative part is then distributed over the surfaces within the zone in some prescribed manner. This, of course, is not a completely realistic model, and it departs from the heat balance principles. However, it is virtually impossible to treat this source in any more detail since the alternative would require knowledge of the placement and surface temperatures of all equipment.
Internal ShortWave Radiation[LINK]
SW Radiation from Lights[LINK]
The short wavelength radiation from lights is distributed over the surfaces in the zone in some prescribed manner.
Transmitted Solar[LINK]
Transmitted solar radiation is also distributed over the surfaces in the zone in a prescribed manner. It would be possible to calculate the actual position of beam solar radiation, but that would involve partial surface irradiation, which is inconsistent with the rest of the zone model that assumes uniform conditions over an entire surface. The current procedures incorporate a set of prescribed distributions. Since the heat balance approach can deal with any distribution function, it is possible to change the distribution function if it seems appropriate.
Convection to Zone Air[LINK]
The convection flux is calculated using the heat transfer coefficients as follows:
q′′conv=hc(Ts−Ta) (89)
The inside convection coefficients (hc) can be calculated using one of many different models. Currently the implementation uses coefficients based on correlations for natural, mixed, and forced convection.
Interior Conduction[LINK]
This contribution to the inside surface heat balance is the wall conduction term, q′′ki shown in Equation (30). This represents the heat transfer to the inside face of the building element. Again, a CTF formulation is used to determine this heat flux.
Interior Convection[LINK]
There are many different modeling options available in EnergyPlus for inside convection coefficients, hc. There are four different settings to direct how EnergyPlus managers select hc models during a simulation. There are numerous individual model equations for hc in EnergyPlus to cover different situations that arise from surface orientations, room airflow conditions, and heat flow direction. Additionally, in many cases multiple researchers have developed competing models for the same situations that differ and there is no way to declare one is better than another. An overall default for the simulation is selected in the “SurfaceConvectionAlgorithm:Inside” object and can be overridden by selecting a different option in a zone description. These models are explained in the following sections. In addition to the correlation choices described below, it is also possible to override the convection coefficients on the inside of any surface by using the “SurfaceProperty:ConvectionCoefficients” object in the input file to set the convection coefficient value on the inside of any surface. The values can be specified directly or with schedules. Specific details are given in the Input Output Reference document.
Adaptive Convection Algorithm[LINK]
BeausoleilMorrison (2000, 2002) developed a methodology for dynamically managing the selection of hc equations called adaptive convection algorithm. The algorithm is used to select among the available hc equations for the one that is most appropriate for a given surface at a given time. As BeausoleilMorrison notes, the adaptive convection algorithm is intended to be expanded and altered to reflect different classification schemes and/or new hc equations. The implementation in EnergyPlus has been modified from the original in the following ways:
An input mechanism is provided (see the
SurfaceConvectionAlgorithm:Inside:AdapativeModelSelections object) so that the user can customize the specific selections of hc equations that are applied for different flow regimes and surface orientations. The changes apply in a general way to the entire model (but can be overridden by setting surface properties).
To avoid requiring additional user input on the position of ZoneHVACtype equipment within a zone, there is no distinction between zones that have convective zone heater equipment located underneath the windows and those that have convective heaters located away from the windows. This applies to the air flow regime associated with convective zone heaters. Using BeausoleilMorrison’s terminology, regimes B1 and B2 are combined into just one B regime.
To avoid requiring additional user input on the position of ZoneHVACtype equipment within a zone, there is no distinction between surfaces that are directly blown on the fan and those that are away from the fan for the air flow regime associated with mechanical circulation from a zone fan (ZoneHVAC type equipment).
The correlation for horizontal free jet developed by Fisher (1995) is not used. Ceiling diffuser models are used for all mechanical circulation from central air system. This decision was made for two reasons: (1) to avoid requiring additional user input on the position of, and momentum generated by, air terminal units, and (2) because Fisher (1995) found that the Coanda effect is so significant that in practice a free horizontal jet is difficult to maintain and mechanicaldriven room airflows generally attach to surfaces and tend to match the flow regime of a ceiling diffuser much more often than a free jet.
EnergyPlus supports arbitrary geometry so surfaces can be tilted with respect to vertical or horizontal. BeausoleilMorrison’s adaptive convection algorithm was originally structured to use hc equations that have no functional dependence on surface tilt angle. However, tilted surfaces do perform differently than vertical or horizontal surface when buoyancy forces are significant. Therefore, the EnergyPlus implementation expands the structure of the algorithm to include additional categories for tilted surfaces. The hc equations developed by Walton (1983) are selected as the defaults for tilted surfaces because they have a functional dependence on tilt angle.
Fohanno and Polidari (2006) produced a new hc equation for vertical walls inside buildings with a simple buoyancy flow regime. They used a theoretical approach based on integral formalism and uniform heat flux (rather than uniform temperature) that covers both laminar and turbulent flow situations. In EnergyPlus, this model is selected as the default in place of the model by Alamdari and Hammond (1983) for vertical walls.
Karadag (2009) produced a new hc equation for ceiling surfaces that are actively chilled. He used computation fluid dynamics and various sized rooms and temperature conditions. In EnergyPlus, this model is selected as the default for surfaces that have active, inceiling cooling (in place of the model by Alamdari and Hammond (1983) for unstable ceilings).
International Standard Organization (ISO) completed Standard 150992003 which includes hc equations for the inside face of windows. EnergyPlus strives to adhere to formal modeling Standards where possible. Therefore the implementation includes a larger structure for the adaptive algorithm that includes additional categories for windows in all flow regimes and ISO 150992003 models are used as the default for windows in natural convection flow regimes. The ISO 15099 model applies to various tilt angles.
Goldstein and Novosalec (2010) produced new hc equations for forced air situations with ceiling slot diffusers along perimeters with significant glazing fractions. They used experiments with fullsized test room. These new equations are selected as the default for windows, ceilings and floors when there is an active central air system.
Interior mass surfaces are assigned the hc equation that would apply (stable or unstable) to a horizontal, upward facing surface for each flow regime.
The algorithm switches between forced, mixed, and natural flow regimes by calculating the Richardson number, Ri = Gr/Re^2, for the zone. Large values of Ri indicate buoyancy dominates, while small values indicate forced flows dominate. To distinguish between opposing Zone unit type equipment (with fans) are assumed to force air up walls, and central air type equipment (with diffusers) are assumed to force air down walls.
The adaptive convection algorithm implemented in EnergyPlus for the inside face has a total of 45 different categories for surfaces and 29 different options for hc equation selections. The following table summarizes the categories and the default assignments for hc equations. The individual hc equations are documented below.
Inside Convection Categories and Assignments
1 


Vertical Walls 
FohannoPolidoriVerticalWall* 




AlamdariHammondVerticalWall 




ASHRAEVerticalWall 
2 


Stable Horizontal 
AlamdariHammondStableHorizontal* 




WaltonStableHorizontalOrTilt 
3 
Simple Buoyancy 
A3 
Unstable Horizontal 
AlamdariHammondUnstableHorizontal* 




WaltonUnstableHorizontalOrTilt 
4 


Stable Tilted 
WaltonStableHorizontalOrTilt* 
5 


Unstable Tilted 
WaltonUnstableHorizontalOrTilt* 
6 


Windows 
ISO15099Windows* 
7 


Vertical Walls 
KhalifaEq3WallAwayFromHeat* 




FohannoPolidoriVerticalWall 




AlamdariHammondVerticalWall 




ASHRAEVerticalWall 
8 


Stable Horizontal 
AlamdariHammondStableHorizontal* 




WaltonStableHorizontalOrTilt 
9 
Infloor Heating or Inceiling Cooling 
A1 
Unstable Horizontal 
KhalifaEq4CeilingAwayFromHeat* 




AlamdariHammondUnstableHorizontal 




WaltonUnstableHorizontalOrTilt 
10 


Heated Floor 
AwbiHattonHeatedFloor* 




WaltonUnstableHorizontalOrTilt 




AlamdariHammondUnstableHorizontal 
11 


Chilled Ceiling 
KaradagChilledCeiling* 




WaltonUnstableHorizontalOrTilt 
12 


Stable Tilted 
WaltonStableHorizontalOrTilt* 
13 


Unstable Tilted 
WaltonUnstableHorizontalOrTilt* 
14 


Windows 
ISO15099Windows* 
15 


Vertical Walls (nonheated) 
KhalifaEq6NonHeatedWalls* 




FohannoPolidoriVerticalWall 




ASHRAEVerticalWall 
16 


Heated Wall 
AwbiHattonHeatedWall* 
17 
Wall Panel Heating 
A2 
Stable Horizontal 
AlamdariHammondStableHorizontal* 




WaltonStableHorizontalOrTilt 
18 


Unstable Horizontal 
KhalifaEq7Ceiling* 




AlamdariHammondUnstableHorizontal 




WaltonUnstableHorizontalOrTilt 
19 


Stable Tilted 
WaltonStableHorizontalOrTilt* 
20 


Unstable Tilted 
WaltonUnstableHorizontalOrTilt* 
21 


Windows 
ISO15099Windows* 
22 


Vertical Walls not near heater 
FohannoPolidoriVerticalWall* 




KhalifaEq6NonHeatedWalls 




KhalifaEq3WallAwayFromHeat 




AlamdariHammondVerticalWall 




ASHRAEVerticalWall 
23 


Vertical Walls near heater 
KhalifaEq5WallNearHeat* 
24 
Convective Zone Heater 
B 
Stable Horizontal 
AlamdariHammondStableHorizontal* 




WaltonStableHorizontalOrTilt 
25 


Unstable Horizontal 
KhalifaEq7Ceiling* 




KhalifaEq4CeilingAwayFromHeat 




WaltonUnstableHorizontalOrTilt 
26 


Stable Tilted 
WaltonStableHorizontalOrTilt* 
27 


Unstable Tilted 
WaltonUnstableHorizontalOrTilt* 
28 


Windows 
ISO15099Windows* 
29 


Walls 
GoldsteinNovoselacCeilingDiffuserWalls* 




FisherPedersenCeilingDiffuserWalls 
30 


Ceiling 
FisherPedersenCeilingDiffuserCeiling* 
31 
Mechanical Central Air Diffuser 
C 
Floor 
GoldsteinNovoselacCeilingDiffuserFloor* 




FisherPedersenCeilingDiffuserFloor 
32 


Windows 
GoldsteinNovoselacCeilingDiffuserWindow* 




ISO15099Windows 
33 


Walls 
KhalifaEq3WallAwayFromHeat * 
34 


Stable Horizontal 
AlamdariHammondStableHorizontal* 




WaltonStableHorizontalOrTilt 
35 
Mechanical Zone Fan Circulation 
D 
Unstable Horizontal 
KhalifaEq4CeilingAwayFromHeat* 




WaltonUnstableHorizontalOrTilt 
36 


Stable Tilted 
WaltonStableHorizontalOrTilt* 
37 


Unstable Tilted 
WaltonUnstableHorizontalOrTilt* 
38 


Windows 
GoldsteinNovoselacCeilingDiffuserWindow* 




ISO15099Windows 
39 


Assisting Flow Walls 
BeausoleilMorrisonMixedAssistedWall* 
40 


Opposing Flow Walls 
BeausoleilMorrisonMixedOpposingWall* 
41 


Stable Floor 
BeausoleilMorrisonMixedStableFloor* 
42 
Mixed 
E 
Unstable Floor 
BeausoleilMorrisonMixedUnstableFloor* 
43 


Stable Ceiling 
BeausoleilMorrisonMixedStableCeiling* 
44 


Unstable Ceiling 
BeausoleilMorrisonMixedUnstableCeiling* 
45 


Windows 
GoldsteinNovoselacCeilingDiffuserWindow* 




ISO15099Windows 
*Indicates the default selection for hc model equation.
Inside Face Surface Classification[LINK]
The adaptive convection algorithm is based on classifying surfaces by flow regime and orientation so that the correct hc equation can be chosen at a particular point in time during the simulation. The classification depends on user input with some aspects processed only once at the beginning and others during each timestep. There are also various parameters or inputs to the hc equations that need static or dynamic processing.
For each surface, it and the zone it is attached to are processed for the following static characteristics.
Characteristic height for convection is taken as the zone height.
Surfaces listed as receiving heat from Zone HVAC equipment with radiative models are considered “near” the heater.
Zones are examined for low temperature radiant systems. The surfaces that contain the active elements are examined and the zone characterized to know if it has infloor heating, inceiling cooling, or inwall heating.
A hydraulic diameter is calculated for horizontal surfaces for the entire zone.
and calculating various parameters needed by hc equations. Selecting flow regime is done in the following manner. For each surface, we examine the zone on the inside face for the following:
HVAC system type
HVAC operating status
HVAC system ACH
The surfaces are evaluated to determine:
Surface classification: Floor, wall, roof, window, types
Tilt angle
Convective stability (sign of ΔT)
The individual hc model equations and their respective references are listed in next by the keyword used to identify them.
ASHRAE Vertical Wall[LINK]
Walton adopted the following equation for natural convection from ASHRAE .
h=1.31ΔT13
This is usually bound at a minimum of .1 in EnergyPlus. This is a component of the TARP overall algorithm described below.
Walton Unstable Horizontal Or Tilt[LINK]
Walton (11983) developed the following equation by fitting curves from various sources.
h=9.482ΔT137.283−cosΣ
Unstable refers to the direction of heat flow and the associated buoyancy relative to the surfaces. Unstable is when the natural tendency is to enhance flow in the sense that rising warmer air, or falling cooler air, is free to move away from the surface. This is usually bound at a minimum of .1 in EnergyPlus. This is a component of the TARP overall algorithm described below.
Walton Stable Horizontal Or Tilt[LINK]
Walton (11983) developed the following equation by fitting curves from various sources.
h=1.810ΔT131.382+cosΣ
Stable refers to the direction of heat flow and the associated buoyancy relative to the surfaces. Stable is when the natural tendency is to retard flow in the sense that rising warmer air, or falling cooler air, is driven against the surface. This is usually bound at a minimum of .1 in EnergyPlus. This is a component of the TARP overall algorithm described below.
Fisher Pedersen Ceiling Diffuser Walls[LINK]
Fisher and Pedersen 1997) developed the following equation from laboratory chamber measurements.
h=1.208+1.012∗ACH0.604
This is a component of the CeilingDiffuser overall algorithm described below.
Fisher Pedersen Ceiling Diffuser Ceiling[LINK]
Fisher and Pedersen (1997) developed the following equation from laboratory chamber measurements.
h=2.234+4.099∗ACH0.503
This is a component of the CeilingDiffuser overall algorithm described below.
Fisher Pedersen Ceiling Diffuser Floor[LINK]
Fisher and Pedersen 1997) developed the following equation from laboratory chamber measurements.
h=3.873+0.082∗ACH0.98
This is a component of the CeilingDiffuser overall algorithm described below.
Alamdari Hammond Stable Horizontal[LINK]
Alamdari and Hammond (1983) developed the following correlation for horizontal surfaces in stable thermal situation.
h=0.6(ΔTD2h)1/155
where,
Dh=4AP , hydraulic diameter of horizontal surface, A is area (m2) and P is the perimeter (m) of the entire zone.
Alamdari Hammond Unstable Horizontal[LINK]
Alamdari and Hammond (1983) developed the following correlation for horizontal surfaces in a buoyant thermal situation.
h=⎧⎪
⎪
⎪
⎪⎨⎪
⎪
⎪
⎪⎩⎡⎢
⎢
⎢⎣1.4(ΔTDh)1/144⎤⎥
⎥
⎥⎦6+[1.63ΔT1/133]6⎫⎪
⎪
⎪
⎪⎬⎪
⎪
⎪
⎪⎭1/166
Alamdari Hammond Vertical Wall[LINK]
Alamdari and Hammond (1983) developed the following correlation for vertical surfaces.
h=⎧⎪
⎪
⎪
⎪⎨⎪
⎪
⎪
⎪⎩⎡⎢
⎢
⎢⎣1.5(ΔTH)1/144⎤⎥
⎥
⎥⎦6+[1.23ΔT1/133]6⎫⎪
⎪
⎪
⎪⎬⎪
⎪
⎪
⎪⎭1/166
where,
H is the characteristic height for the surface. In EnergyPlus this is the zone’s ceiling height (which could be larger than the height of an individual surface when wall are subdivided into more than one surface).
Khalifa Eq3 Wall Away From Heat[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 3” in original reference, is for convectively heated zones and applies to the inside surfaces of walls away from the heat source:
h=2.07ΔT0.23
Khalifa Eq4 Ceiling Away From Heat[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 4” in original reference, is for convectively heated zones and applies to the inside surfaces of ceilings away from the heat source:
h=2.72ΔT0.13
Khalifa Eq5 Wall Near Heat[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 5” in original reference, is for convectively heated zones and applies to the inside surfaces of walls near the heat source:
h=1.98ΔT0.32
Khalifa Eq6 Non Heated Walls[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 6” in original reference, is for heated zones and applies to the inside surfaces of walls that are not heated:
h=2.30ΔT0.24
Khalifa Eq7 Ceiling[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 7” in original reference, is for heated zones and applies to the inside surfaces of ceilings:
h=3.10ΔT0.17
Awbi Hatton Heated Floor[LINK]
Awbi and Hatton (1999) conducted laboratory measurements using environmental chambers and developed the following correlation for floor surfaces that are being actively heated.
h=2.175ΔT0.308D0.076h
where,
Dh=4AP, hydraulic diameter of horizontal surface, A is area (m2) and P is the perimeter (m) of the entire zone (all of the adjacent floor surfaces if more than one in the zone).
Awbi Hatton Heated Wall[LINK]
Awbi and Hatton (1999) developed the following correlation for wall surfaces that are being actively heated.
h=1.823ΔT0.293D0.076h
where,
Dh=4AP, hydraulic diameter of wall surface, A is area (m2) and P is the perimeter (m) of the entire wall (all of the adjacent wall surfaces if more than one along the wall).
Beausoleil Morrison Mixed Assisted Wall[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for walls where the flow driving forces from mechanical forces are augmented by the driving forces from buoyancy.
h=⎛⎜
⎜⎝⎧⎪⎨⎪⎩⎡⎢⎣1.5(ΔTH)1/4⎤⎥⎦6+[1.23ΔT2]1/6⎫⎪⎬⎪⎭1/2+{[Tsurf−TSATΔT][−0.199+0.190⋅ACH0.8]}3⎞⎟
⎟⎠1/3
where,
TSAT is the supply air temperature at the diffuser.
Here the reference temperature is the zone air temperature rather than the diffuser supply air temperature.
Beausoleil Morrison Mixed Opposing Wall[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for walls where the flow driving forces from mechanical forces are opposed by the driving forces from buoyancy.
h=max⎧⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪⎨⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪⎩⎛⎝{[1.5(ΔTH)1/4]6+[1.23ΔT2]1/6}1/2−{[Tsurf−TSATΔT]⋅[−0.199+0.190⋅ACH0.8]}3⎞⎠1/30.8⋅{[1.5(ΔTH)1/4]6+[1.23ΔT2]}1/60.8⋅{[Tsurf−TSATΔT]⋅[−0.199+0.190⋅ACH0.8]}
Beausoleil Morrison Mixed Stable Floor[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for floors where the flow driving forces include both mechanical forces and thermally stable buoyancy.
h=⎛⎜
⎜⎝⎧⎪
⎪⎨⎪
⎪⎩0.6⋅(ΔTDH)1/155⎫⎪
⎪⎬⎪
⎪⎭3+{[Tsurf−TSATΔT]⋅[0.159+0.116ACH0.8]}3⎞⎟
⎟⎠1/133
Beausoleil Morrison Mixed Unstable Floor[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for floors where the flow driving forces include both mechanical forces and thermally unstable buoyancy.
h=⎛⎜
⎜
⎜
⎜
⎜
⎜⎝⎧⎪
⎪
⎪
⎪⎨⎪
⎪
⎪
⎪⎩⎡⎢
⎢
⎢⎣1.4(ΔTDh)1/144⎤⎥
⎥
⎥⎦6+[1.63ΔT1/133]6⎫⎪
⎪
⎪
⎪⎬⎪
⎪
⎪
⎪⎭3/366+{[Tsurf−TSATΔT]⋅[0.159+0.116ACH0.8]}3⎞⎟
⎟
⎟
⎟
⎟
⎟⎠1/133
Beausoleil Morrison Mixed Stable Ceiling[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for ceilings where the flow driving forces include both mechanical forces and thermally Stable buoyancy.
h=⎛⎜
⎜⎝⎧⎪
⎪⎨⎪
⎪⎩0.6⋅(ΔTDH)1/155⎫⎪
⎪⎬⎪
⎪⎭3+{[Tsurf−TSATΔT]⋅[−0.166+0.484ACH0.8]}3⎞⎟
⎟⎠1/133
Beausoleil Morrison Mixed Unstable Ceiling[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for ceilings where the flow driving forces include both mechanical forces and thermally unstable buoyancy.
h=⎛⎜
⎜⎝⎛⎜⎝⎡⎢⎣1.4(ΔTDh)1/4⎤⎥⎦6+[1.63ΔT1/3]6⎞⎟⎠3/6+([Tsurf−TSATΔT]⋅[−0.166+0.484⋅ACH0.8])3⎞⎟
⎟⎠1/3
Fohanno Polidori Vertical Wall[LINK]
Fohanno and Polidori (2006) developed the following equation for hc for vertical walls under simple buoyancy flow conditions.
h=⎧⎪
⎪⎨⎪
⎪⎩1.332(ΔTH)1/144,Ra∗H≤6.3×1091.235e(0.0467H)ΔT0.316,Ra∗H>6.3×109
where,
Ra∗H=gβfq′′cH4kfν2fPrf
Karadag Chilled Ceiling[LINK]
Karadag (2009) used numerical methods to develop the following equation for hc for chilledceiling surfaces.
h=3.1ΔT0.22
ISO 15099 Windows[LINK]
ISO Standard 150992003 includes the equations for hc for roomside of windows and surfaces with any tilt angle and heat flow direction. The ISO 15099 correlation is for still room air angle and is determined in terms of the Nusselt number, Nu , where
hi=Nu(λH)
where,
λ is the thermal conductivity of air, and
H is the height of the window.
The Rayleigh number based on height, RaH , is calculated using,
RaH=ρ2H3gcp∣∣Tsurf,i−Tair∣∣Tm,fμλ
where,
ρ is the density of air
g is the acceleration due to gravity,
cp is the specific heat of air,
μ is the dynamic viscosity of air, and
Tm,f is the mean film temperature in Kelvin given by,
Tm,f=Tair+14(Tsurf,i−Tair)
There are four cases for the Nusselt correlation that vary by the tilt angle in degrees, γ , and are based on heating conditions. For cooling conditions (where Tsurf,i>Tair ) the tilt angle is complemented so that γ=180−γ
Case A. 0o≤γ<15o
Nu=0.13Ra1/133H
Case B. 15o≤γ≤90o
Racv=2.5×105(e0.72γsinλ)1/155
Nu=0.56(RaHsinγ)1/144;forRaH≤RaCV
Nu=0.13(Ra1/133H−Ra1/133CV)+0.56(RaCVsinγ)1/144;RaH>RaCV
Case C. 90o<γ≤179o
Nu=0.56(RaHsinγ)1/144;105≤RaHsinγ<1011
Case D. 179o<γ≤180o
Nu=0.58Ra1/155H;RaH≤1011
The material properties are evaluated at the mean film temperature. Standard EnergyPlus psychrometric functions are used for ρ and cp . Thermal conductivity is calculated using,
λ=2.873×10−3+7.76×10−8Tm,f .
Kinematic viscosity is calculated using,
μ=3.723×10−6+4.94×10−8Tm,f .
This correlation depends on the surface temperature of the roomside glazing surface and is therefore included inside the window heat balance iteration loop.
Goldstein Novoselac Ceiling Diffuser Window[LINK]
Goldstein and Novoselac (2010) used laboratory chamber measurements to develop convection correlations for perimeter zones with highly glazed spaces served by overhead slotdiffuserbased air systems. The following are for bare windows in such spaces.
For WWR<50% with window in upper part of wall:
h=0.117(˙VL)0.8
For WWR<50% with window in lower part of wall
h=0.093(˙VL)0.8
For WWR > 50%
h=0.103(˙VL)0.8
Where,
WWR is the window to wall ratio.
L is the length of exterior wall with glazing in the zone.
˙V is the air system flow rate in m3/s.
Goldstein Novoselac Ceiling Diffuser Walls[LINK]
Goldstein and Novoselac (2010) used laboratory chamber measurements to develop convection correlations for perimeter zones with highly glazed spaces served by overhead slotdiffuserbased air systems. The following are for exterior walls in such spaces.
For walls located below a window
h=0.063(˙VL)0.8
For walls located above a window
h=0.093(˙VL)0.8
Goldstein Novoselac Ceiling Diffuser Floor[LINK]
Goldstein and Novoselac (2010) used laboratory chamber measurements to develop convection correlations for perimeter zones with highly glazed spaces served by overhead slotdiffuserbased air systems. The following is for floors in such spaces.
h=0.048(˙VL)0.8
Separate from the above model structure, there are also other comprehensive algorithm structures which are described below.
TARP Algorithm[LINK]
The comprehensive natural convection model, accessed using the keyword “TARP,” correlates the convective heat transfer coefficient to the surface orientation and the difference between the surface and zone air temperatures (where DT = Surface Temperature  Air Temperature). The algorithm is taken directly from Walton (1983). Walton derived his algorithm from ASHRAE literature which can now be found for example in the ASHRAE Handbook (HoF 2001), Table 5 on p. 3.12, which gives equations for natural convection heat transfer coefficients in the turbulent range for large, vertical plates and for large, horizontal plates facing upward when heated (or downward when cooled). A note in the text also gives an approximation for large, horizontal plates facing downward when heated (or upward when cooled) recommending that it should be half of the facing upward value. Walton adds a curve fit as a function of the cosine of the tilt angle to provide intermediate values between vertical and horizontal. The curve fit values at the extremes match the ASHRAE values very well.
For no temperature difference OR a vertical surface the following correlation is used:
h=1.31ΔT13 (90)
For (ΔT < 0.0 AND an upward facing surface) OR (ΔT > 0.0 AND an downward facing surface) an enhanced convection correlation is used:
h=9.482ΔT137.283−cosΣ (91)
where Σ is the surface tilt angle.
For (ΔT > 0.0 AND an upward facing surface) OR (ΔT < 0.0 AND an downward facing surface) a reduced convection correlation is used:
h=1.810ΔT131.382+cosΣ (92)
where Σ is the surface tilt angle.
Simple Natural Convection Algorithm[LINK]
The simple convection model uses constant coefficients for different heat transfer configurations, using the same criteria as the detailed model to determine reduced and enhanced convection. The coefficients are also taken directly from Walton (1983). Walton derived his coefficients from the surface conductances for e = 0.90 found in the ASHRAE Handbook (1985) in Table on p. 23.2. The radiative heat transfer component was estimated at 1.02 * 0.9 = 0.918 BTU/hft2F and then subtracted off. Finally the coefficients were converted to SI units to yield the values below.
For a vertical surface:
h=3.076
For a horizontal surface with reduced convection:
h=0.948
For a horizontal surface with enhanced convection:
h=4.040
For a tilted surface with reduced convection:
h=2.281
For a tilted surface with enhanced convection:
h=3.870
Ceiling Diffuser Algorithm[LINK]
The ceiling diffuser algorithm is based on empirical correlations developed by Fisher and Pedersen (1997). The correlation was reformulated to use the room outlet temperature as the reference temperature. The correlations are shown below.
For Floors:
h=3.873+0.082∗ACH0.98 (93)
The correlation for floors is illustrated in the following figure:
For ceilings:
h=2.234+4.099∗ACH0.503 (94)
The correlation for ceilings is illustrated in the following figure:
For Walls:
h=1.208+1.012∗ACH0.604 (95)
The correlation for walls is illustrated in the following figure:
Trombe Wall Algorithm[LINK]
The Trombe wall algorithm is used to model convection in a “Trombe wall zone”, i.e. the air space between the storage wall surface and the exterior glazing. (See the later sections on Passive and Active Trombe Walls below for more information about Trombe walls.) The algorithm is identical to the convection model (based on ISO 15099) used in Window5 for convection between glazing layers in multipane window systems. The use of the algorithm for modeling an unvented Trombe wall has been validated against experimental data by Ellis (2003).
This algorithm gives the convection coefficients for air in a narrow vertical cavity that is sealed and not ventilated. This applies both to the air gap in between panes of a window or to the air gap between the Trombe wall glazing and the inner surface (often a selective surface). These convection coefficients are really the only difference between a normal zone and a Trombe zone. The rest of the zone heat balance is the same, e.g., transmitted solar, longwave radiation between surfaces, etc.
For a vertical cavity, the correlation from ISO 15099 is:
NU1=0.0673838Ra13 for 5E4 < Ra < 1E6
NU1=0.028154Ra0.4134 for 1E4 < Ra < 5E4
NU1=1+1.7596678E−10Ra2.2984755 for Ra < = 1E4
NU2=0.242({{\rm{Ra}}}/RaA{\rm{A}})0.272
NU=MAX(NU1,NU2)
where
Nu = Nusselt number
Ra = Rayleigh number
A = aspect ratio of cavity
This is then used in EnergyPlus as follows:
Net convection coefficient from glazing to wall is:
hnet=k({{\rm{NU}}}/NUL{\rm{L}})
where
k = conductivity of air
L = air gap thickness
Convection coefficient applied to each wall separately and actually used in the zone heat balance is:
hc=2hnet
Alamdari, F. and G.P. Hammond. 1983. Improved data correlations for buoyancydriven convection in rooms. Building Services Engineering Research & Technology. Vol. 4, No. 3.
ASHRAE. 1985. 1985 ASHRAE Handbook – Fundamentals, Atlanta: American Society of Heating, Refrigerating, and AirConditioning Engineers, Inc.
ASHRAE. 2001. 2001 ASHRAE Handbook – Fundamentals, Atlanta: American Society of Heating, Refrigerating, and AirConditioning Engineers, Inc.
Awbi, H.B. and A. Hatton. 1999. Natural convection from heated room surfaces. Energy and Buildings 30 (1999) 233244.
BeausoleilMorrison, I. 2000. The adaptive coupling of heat and air flow modeling within dynamic wholebuilding simulations. PhD. Thesis. University of Strathclyde, Glasgow, UK.
Ellis, Peter G. 2003. Development and Validation of the Unvented Trombe Wall Model in EnergyPlus. Master’s Thesis, University of Illinois at UrbanaChampaign.
Fisher, D.E. and C.O. Pedersen. 1997. “Convective Heat Transfer in Building Energy and Thermal Load Calculations”, ASHRAE Transactions, Vol. 103, Pt. 2.
Fohanno, S., and G. Polidori. 2006. Modelling of natural convective heat transfer at an internal surface. Energy and Buildings 38 (2006) 548  553
Goldstein, K. and A. Novoselac. 2010. Convective Heat Transfer in Rooms With Ceiling Slot Diffusers (RP1416). HVAC&R Research Journal TBD
Karadag, R. 2009. New approach relevant to total heat transfer coefficient including the effect of radiation and convection at the ceiling in a cooled ceiling room. Applied Thermal Engineering 29 (2009) 15611565
Khalifa AJN. 1989. Heat transfer processes in buildings. Ph.D. Thesis, University of Wales College of Cardiff, Cardiff, UK.
ISO. 2003. ISO 15099:2003. Thermal performance of windows, doors, and shading devices – Detailed calculations. International Organization for Standardization.
Walton, G. N. 1983. Thermal Analysis Research Program Reference Manual. NBSSIR 832655. National Bureau of Standards (now NIST). This is documentation for “TARP.”
Inside Heat Balance[LINK]
The heart of the heat balance method is the internal heat balance involving the inside faces of the zone surfaces. This heat balance is generally modeled with four coupled heat transfer components: 1) conduction through the building element, 2) convection to the air, 3) short wave radiation absorption and reflectance and 4) longwave radiant interchange. The incident short wave radiation is from the solar radiation entering the zone through windows and emittance from internal sources such as lights. The longwave radiation interchange includes the absorption and emittance of low temperature radiation sources, such as all other zone surfaces, equipment, and people.
The heat balance on the inside face can be written as follows:
q′′LWX+q′′SW+q′′LWS+q′′ki+q′′sol+q′′conv=0
where:
q′′LWX = Net longwave radiant exchange flux between zone surfaces.
q′′SW = Net short wave radiation flux to surface from lights.
q′′LWS = Longwave radiation flux from equipment in zone.
q′′ki = Conduction flux through the wall.
q′′sol = Transmitted solar radiation flux absorbed at surface.
q′′conv = Convective heat flux to zone air.
Each of these heat balance components is introduced briefly below.
Inside Heat Balance Control Volume Diagram
Internal LongWave Radiation Exchange[LINK]
LW Radiation Exchange Among Zone Surfaces[LINK]
There are two limiting cases for internal LW radiation exchange that are easily modeled:
The zone air is completely transparent to LW radiation.
The zone air completely absorbs LW radiation from the surfaces within the zone.
The limiting case of completely absorbing air has been used for load calculations and also in some energy analysis calculations. This model is attractive because it can be formulated simply using a combined radiation and convection heat transfer coefficient from each surface to the zone air. However, it oversimplifies the zone surface exchange problem, and as a result, the heat balance formulation in EnergyPlus treats air as completely transparent. This means that it does not participate in the LW radiation exchange among the surfaces in the zone. The model, which considers room air to be completely transparent, is reasonable physically because of the low water vapor concentrations and the short mean path lengths. It also permits separating the radiant and convective parts of the heat transfer at the surface, which is an important attribute of the heat balance method.
EnergyPlus uses a grey interchange model for the longwave radiation among zone surfaces. This model is based on the “ScriptF” concept developed by Hottel (Hottel and Sarofim, Radiative Transfer, Chapter 3, McGraw Hill, 1967). This procedure relies on a matrix of exchange coefficients between pairs of surfaces that include all exchange paths between the surfaces. In other words all reflections, absorptions and reemissions from other surfaces in the enclosure are included in the exchange coefficient, which is called ScriptF. The major assumptions are that all surface radiation properties are grey and all radiation is diffuse. Both assumptions are reasonable for building zone interchange.
The ScriptF coefficients are developed by starting with the traditional direct radiation view factors. In the case of building rooms and zones, there are several complicating factors in finding the direct view factors—the main one being that the location of surfaces such as thermal mass representing furniture and partitions are not known. The other limitation is that the exact calculation of direct view factors is computationally very intensive even if the positions of all surfaces are known. Accordingly, EnergyPlus uses a procedure to approximate the direct view factors. The procedure has two steps:
1) Determine the total area of other surfaces “seen” by a surface.
2) Approximate the direct view factor from surface 1 to surface 2 as the ratio of the area of surface 2 to the total area “seen” by surface 1.
The determination of the “seen” area has several constraints:
No surface sees itself.
All surfaces see thermal mass surfaces.
No surface facing within 10 degrees of another surface is seen by the other surface.
All surfaces see roofs, floors and ceilings (subject to the preceding facing direction constraint).
Because the approximate view factors may not satisfy the basic requirements of reciprocity (two surfaces should exchange equal amounts of heat in each direction), and completeness (every surface should have a direct view factor sum of 1.0), EnergyPlus does a view factor fix operation before they are used in the ScriptF determination. Normally both of the requirements are satisfied, but in some special situations they are not, and special rules are applied.
If a user includes less than four surfaces in a zone, the enforcement of “completeness” is not strictly possible because it is not possible to have a completed enclosure with less than four surfaces. In this situation, reciprocity is enforced initially and then view factors from each surface are summed. When dealing with this few surfaces, the result of using approximate view factors is that once reciprocity is enforced that the largest surface in the zone will end up having view factors that sum up to greater than unity. This is not physically possible. So, to correct this and maintain completeness for one surface, the entire view factor matrix is divided by the largest summation of view factors for any surface in the zone. By dividing the entire view factor matrix by the same value, the reciprocity established by the standard algorithm is maintained and completeness for one surface is established. All of the other surfaces will not have view factors that sum up to unity and thus will not satisfy completeness. However, these other surfaces will not have view factors that sum up to greater than unity thus avoiding a physically impossible situation.
If the area of one surface in a zone is greater than the sum of the areas of all other surfaces, reciprocity only is enforced, but sometimes, for very large surfaces, that enforcement becomes impossible, and the view factors are modified so that only the large surface is seen by very small surfaces.
Warning messages are produced for both of these cases, and the results should be examined very carefully to ascertain that they are reasonable. The suggested action for the second case (the extralarge surface) is to divide the large surface into several smaller surfaces; then the enclosure will be treated as normal.
Once the ScriptF coefficients are determined, the longwave radiant exchange is calculated for each surface using:
qi,j=AiFi,j(T4i−T4j)
where Fi,j is the ScriptF between surfaces i and j.
Thermal Mass and Furniture[LINK]
Furniture in a zone has the effect of increasing the amount of surface area that can participate in the radiation and convection heat exchanges. It also adds participating thermal mass to the zone. These two changes both affect the response to temperature changes in the zone and also affect the heat extraction characteristics.
The proper modeling of furniture is an area that needs further research, but the heat balance formulation allows the effect to be modeled in a realistic manner by including the furniture surface area and thermal mass in the heat exchange process.
LW Radiation From Internal Sources[LINK]
The traditional model for this source is to define a radiative/convective split for the heat introduced into a zone from equipment. The radiative part is then distributed over the surfaces within the zone in some prescribed manner. This, of course, is not a completely realistic model, and it departs from the heat balance principles. However, it is virtually impossible to treat this source in any more detail since the alternative would require knowledge of the placement and surface temperatures of all equipment.
Internal ShortWave Radiation[LINK]
SW Radiation from Lights[LINK]
The short wavelength radiation from lights is distributed over the surfaces in the zone in some prescribed manner.
Transmitted Solar[LINK]
Transmitted solar radiation is also distributed over the surfaces in the zone in a prescribed manner. It would be possible to calculate the actual position of beam solar radiation, but that would involve partial surface irradiation, which is inconsistent with the rest of the zone model that assumes uniform conditions over an entire surface. The current procedures incorporate a set of prescribed distributions. Since the heat balance approach can deal with any distribution function, it is possible to change the distribution function if it seems appropriate.
Convection to Zone Air[LINK]
The convection flux is calculated using the heat transfer coefficients as follows:
q′′conv=hc(Ts−Ta) (89)
The inside convection coefficients (hc) can be calculated using one of many different models. Currently the implementation uses coefficients based on correlations for natural, mixed, and forced convection.
Interior Conduction[LINK]
This contribution to the inside surface heat balance is the wall conduction term, q′′ki shown in Equation (30). This represents the heat transfer to the inside face of the building element. Again, a CTF formulation is used to determine this heat flux.
Interior Convection[LINK]
There are many different modeling options available in EnergyPlus for inside convection coefficients, hc. There are four different settings to direct how EnergyPlus managers select hc models during a simulation. There are numerous individual model equations for hc in EnergyPlus to cover different situations that arise from surface orientations, room airflow conditions, and heat flow direction. Additionally, in many cases multiple researchers have developed competing models for the same situations that differ and there is no way to declare one is better than another. An overall default for the simulation is selected in the “SurfaceConvectionAlgorithm:Inside” object and can be overridden by selecting a different option in a zone description. These models are explained in the following sections. In addition to the correlation choices described below, it is also possible to override the convection coefficients on the inside of any surface by using the “SurfaceProperty:ConvectionCoefficients” object in the input file to set the convection coefficient value on the inside of any surface. The values can be specified directly or with schedules. Specific details are given in the Input Output Reference document.
Adaptive Convection Algorithm[LINK]
BeausoleilMorrison (2000, 2002) developed a methodology for dynamically managing the selection of hc equations called adaptive convection algorithm. The algorithm is used to select among the available hc equations for the one that is most appropriate for a given surface at a given time. As BeausoleilMorrison notes, the adaptive convection algorithm is intended to be expanded and altered to reflect different classification schemes and/or new hc equations. The implementation in EnergyPlus has been modified from the original in the following ways:
An input mechanism is provided (see the
SurfaceConvectionAlgorithm:Inside:AdapativeModelSelections object) so that the user can customize the specific selections of hc equations that are applied for different flow regimes and surface orientations. The changes apply in a general way to the entire model (but can be overridden by setting surface properties).
To avoid requiring additional user input on the position of ZoneHVACtype equipment within a zone, there is no distinction between zones that have convective zone heater equipment located underneath the windows and those that have convective heaters located away from the windows. This applies to the air flow regime associated with convective zone heaters. Using BeausoleilMorrison’s terminology, regimes B1 and B2 are combined into just one B regime.
To avoid requiring additional user input on the position of ZoneHVACtype equipment within a zone, there is no distinction between surfaces that are directly blown on the fan and those that are away from the fan for the air flow regime associated with mechanical circulation from a zone fan (ZoneHVAC type equipment).
The correlation for horizontal free jet developed by Fisher (1995) is not used. Ceiling diffuser models are used for all mechanical circulation from central air system. This decision was made for two reasons: (1) to avoid requiring additional user input on the position of, and momentum generated by, air terminal units, and (2) because Fisher (1995) found that the Coanda effect is so significant that in practice a free horizontal jet is difficult to maintain and mechanicaldriven room airflows generally attach to surfaces and tend to match the flow regime of a ceiling diffuser much more often than a free jet.
EnergyPlus supports arbitrary geometry so surfaces can be tilted with respect to vertical or horizontal. BeausoleilMorrison’s adaptive convection algorithm was originally structured to use hc equations that have no functional dependence on surface tilt angle. However, tilted surfaces do perform differently than vertical or horizontal surface when buoyancy forces are significant. Therefore, the EnergyPlus implementation expands the structure of the algorithm to include additional categories for tilted surfaces. The hc equations developed by Walton (1983) are selected as the defaults for tilted surfaces because they have a functional dependence on tilt angle.
Fohanno and Polidari (2006) produced a new hc equation for vertical walls inside buildings with a simple buoyancy flow regime. They used a theoretical approach based on integral formalism and uniform heat flux (rather than uniform temperature) that covers both laminar and turbulent flow situations. In EnergyPlus, this model is selected as the default in place of the model by Alamdari and Hammond (1983) for vertical walls.
Karadag (2009) produced a new hc equation for ceiling surfaces that are actively chilled. He used computation fluid dynamics and various sized rooms and temperature conditions. In EnergyPlus, this model is selected as the default for surfaces that have active, inceiling cooling (in place of the model by Alamdari and Hammond (1983) for unstable ceilings).
International Standard Organization (ISO) completed Standard 150992003 which includes hc equations for the inside face of windows. EnergyPlus strives to adhere to formal modeling Standards where possible. Therefore the implementation includes a larger structure for the adaptive algorithm that includes additional categories for windows in all flow regimes and ISO 150992003 models are used as the default for windows in natural convection flow regimes. The ISO 15099 model applies to various tilt angles.
Goldstein and Novosalec (2010) produced new hc equations for forced air situations with ceiling slot diffusers along perimeters with significant glazing fractions. They used experiments with fullsized test room. These new equations are selected as the default for windows, ceilings and floors when there is an active central air system.
Interior mass surfaces are assigned the hc equation that would apply (stable or unstable) to a horizontal, upward facing surface for each flow regime.
The algorithm switches between forced, mixed, and natural flow regimes by calculating the Richardson number, Ri = Gr/Re^2, for the zone. Large values of Ri indicate buoyancy dominates, while small values indicate forced flows dominate. To distinguish between opposing Zone unit type equipment (with fans) are assumed to force air up walls, and central air type equipment (with diffusers) are assumed to force air down walls.
The adaptive convection algorithm implemented in EnergyPlus for the inside face has a total of 45 different categories for surfaces and 29 different options for hc equation selections. The following table summarizes the categories and the default assignments for hc equations. The individual hc equations are documented below.
*Indicates the default selection for hc model equation.
Inside Face Surface Classification[LINK]
The adaptive convection algorithm is based on classifying surfaces by flow regime and orientation so that the correct hc equation can be chosen at a particular point in time during the simulation. The classification depends on user input with some aspects processed only once at the beginning and others during each timestep. There are also various parameters or inputs to the hc equations that need static or dynamic processing.
For each surface, it and the zone it is attached to are processed for the following static characteristics.
Characteristic height for convection is taken as the zone height.
Surfaces listed as receiving heat from Zone HVAC equipment with radiative models are considered “near” the heater.
Zones are examined for low temperature radiant systems. The surfaces that contain the active elements are examined and the zone characterized to know if it has infloor heating, inceiling cooling, or inwall heating.
A hydraulic diameter is calculated for horizontal surfaces for the entire zone.
and calculating various parameters needed by hc equations. Selecting flow regime is done in the following manner. For each surface, we examine the zone on the inside face for the following:
HVAC system type
HVAC operating status
HVAC system ACH
The surfaces are evaluated to determine:
Surface classification: Floor, wall, roof, window, types
Tilt angle
Convective stability (sign of ΔT)
The individual hc model equations and their respective references are listed in next by the keyword used to identify them.
ASHRAE Vertical Wall[LINK]
Walton adopted the following equation for natural convection from ASHRAE .
h=1.31ΔT13
This is usually bound at a minimum of .1 in EnergyPlus. This is a component of the TARP overall algorithm described below.
Walton Unstable Horizontal Or Tilt[LINK]
Walton (11983) developed the following equation by fitting curves from various sources.
h=9.482ΔT137.283−cosΣ
Unstable refers to the direction of heat flow and the associated buoyancy relative to the surfaces. Unstable is when the natural tendency is to enhance flow in the sense that rising warmer air, or falling cooler air, is free to move away from the surface. This is usually bound at a minimum of .1 in EnergyPlus. This is a component of the TARP overall algorithm described below.
Walton Stable Horizontal Or Tilt[LINK]
Walton (11983) developed the following equation by fitting curves from various sources.
h=1.810ΔT131.382+cosΣ
Stable refers to the direction of heat flow and the associated buoyancy relative to the surfaces. Stable is when the natural tendency is to retard flow in the sense that rising warmer air, or falling cooler air, is driven against the surface. This is usually bound at a minimum of .1 in EnergyPlus. This is a component of the TARP overall algorithm described below.
Fisher Pedersen Ceiling Diffuser Walls[LINK]
Fisher and Pedersen 1997) developed the following equation from laboratory chamber measurements.
h=1.208+1.012∗ACH0.604
This is a component of the CeilingDiffuser overall algorithm described below.
Fisher Pedersen Ceiling Diffuser Ceiling[LINK]
Fisher and Pedersen (1997) developed the following equation from laboratory chamber measurements.
h=2.234+4.099∗ACH0.503
This is a component of the CeilingDiffuser overall algorithm described below.
Fisher Pedersen Ceiling Diffuser Floor[LINK]
Fisher and Pedersen 1997) developed the following equation from laboratory chamber measurements.
h=3.873+0.082∗ACH0.98
This is a component of the CeilingDiffuser overall algorithm described below.
Alamdari Hammond Stable Horizontal[LINK]
Alamdari and Hammond (1983) developed the following correlation for horizontal surfaces in stable thermal situation.
h=0.6(ΔTD2h)1/155
where,
Dh=4AP , hydraulic diameter of horizontal surface, A is area (m2) and P is the perimeter (m) of the entire zone.
Alamdari Hammond Unstable Horizontal[LINK]
Alamdari and Hammond (1983) developed the following correlation for horizontal surfaces in a buoyant thermal situation.
h=⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩⎡⎢ ⎢ ⎢⎣1.4(ΔTDh)1/144⎤⎥ ⎥ ⎥⎦6+[1.63ΔT1/133]6⎫⎪ ⎪ ⎪ ⎪⎬⎪ ⎪ ⎪ ⎪⎭1/166
Alamdari Hammond Vertical Wall[LINK]
Alamdari and Hammond (1983) developed the following correlation for vertical surfaces.
h=⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩⎡⎢ ⎢ ⎢⎣1.5(ΔTH)1/144⎤⎥ ⎥ ⎥⎦6+[1.23ΔT1/133]6⎫⎪ ⎪ ⎪ ⎪⎬⎪ ⎪ ⎪ ⎪⎭1/166
where,
H is the characteristic height for the surface. In EnergyPlus this is the zone’s ceiling height (which could be larger than the height of an individual surface when wall are subdivided into more than one surface).
Khalifa Eq3 Wall Away From Heat[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 3” in original reference, is for convectively heated zones and applies to the inside surfaces of walls away from the heat source:
h=2.07ΔT0.23
Khalifa Eq4 Ceiling Away From Heat[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 4” in original reference, is for convectively heated zones and applies to the inside surfaces of ceilings away from the heat source:
h=2.72ΔT0.13
Khalifa Eq5 Wall Near Heat[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 5” in original reference, is for convectively heated zones and applies to the inside surfaces of walls near the heat source:
h=1.98ΔT0.32
Khalifa Eq6 Non Heated Walls[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 6” in original reference, is for heated zones and applies to the inside surfaces of walls that are not heated:
h=2.30ΔT0.24
Khalifa Eq7 Ceiling[LINK]
Khalifa (1989) conducted experiments with test chambers and developed correlations for certain types of surfaces. One of them, identified as “Equation 7” in original reference, is for heated zones and applies to the inside surfaces of ceilings:
h=3.10ΔT0.17
Awbi Hatton Heated Floor[LINK]
Awbi and Hatton (1999) conducted laboratory measurements using environmental chambers and developed the following correlation for floor surfaces that are being actively heated.
h=2.175ΔT0.308D0.076h
where,
Dh=4AP, hydraulic diameter of horizontal surface, A is area (m2) and P is the perimeter (m) of the entire zone (all of the adjacent floor surfaces if more than one in the zone).
Awbi Hatton Heated Wall[LINK]
Awbi and Hatton (1999) developed the following correlation for wall surfaces that are being actively heated.
h=1.823ΔT0.293D0.076h
where,
Dh=4AP, hydraulic diameter of wall surface, A is area (m2) and P is the perimeter (m) of the entire wall (all of the adjacent wall surfaces if more than one along the wall).
Beausoleil Morrison Mixed Assisted Wall[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for walls where the flow driving forces from mechanical forces are augmented by the driving forces from buoyancy.
h=⎛⎜ ⎜⎝⎧⎪⎨⎪⎩⎡⎢⎣1.5(ΔTH)1/4⎤⎥⎦6+[1.23ΔT2]1/6⎫⎪⎬⎪⎭1/2+{[Tsurf−TSATΔT][−0.199+0.190⋅ACH0.8]}3⎞⎟ ⎟⎠1/3
where,
TSAT is the supply air temperature at the diffuser.
Here the reference temperature is the zone air temperature rather than the diffuser supply air temperature.
Beausoleil Morrison Mixed Opposing Wall[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for walls where the flow driving forces from mechanical forces are opposed by the driving forces from buoyancy.
h=max⎧⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎩⎛⎝{[1.5(ΔTH)1/4]6+[1.23ΔT2]1/6}1/2−{[Tsurf−TSATΔT]⋅[−0.199+0.190⋅ACH0.8]}3⎞⎠1/30.8⋅{[1.5(ΔTH)1/4]6+[1.23ΔT2]}1/60.8⋅{[Tsurf−TSATΔT]⋅[−0.199+0.190⋅ACH0.8]}
Beausoleil Morrison Mixed Stable Floor[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for floors where the flow driving forces include both mechanical forces and thermally stable buoyancy.
h=⎛⎜ ⎜⎝⎧⎪ ⎪⎨⎪ ⎪⎩0.6⋅(ΔTDH)1/155⎫⎪ ⎪⎬⎪ ⎪⎭3+{[Tsurf−TSATΔT]⋅[0.159+0.116ACH0.8]}3⎞⎟ ⎟⎠1/133
Beausoleil Morrison Mixed Unstable Floor[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for floors where the flow driving forces include both mechanical forces and thermally unstable buoyancy.
h=⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩⎡⎢ ⎢ ⎢⎣1.4(ΔTDh)1/144⎤⎥ ⎥ ⎥⎦6+[1.63ΔT1/133]6⎫⎪ ⎪ ⎪ ⎪⎬⎪ ⎪ ⎪ ⎪⎭3/366+{[Tsurf−TSATΔT]⋅[0.159+0.116ACH0.8]}3⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠1/133
Beausoleil Morrison Mixed Stable Ceiling[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for ceilings where the flow driving forces include both mechanical forces and thermally Stable buoyancy.
h=⎛⎜ ⎜⎝⎧⎪ ⎪⎨⎪ ⎪⎩0.6⋅(ΔTDH)1/155⎫⎪ ⎪⎬⎪ ⎪⎭3+{[Tsurf−TSATΔT]⋅[−0.166+0.484ACH0.8]}3⎞⎟ ⎟⎠1/133
Beausoleil Morrison Mixed Unstable Ceiling[LINK]
BeausoleilMorrison (2000) used blending techniques to combine correlations originally developed by Alamdari and Hammond (1983) and Fisher and Pedersen (1997) to create the following correlation is for ceilings where the flow driving forces include both mechanical forces and thermally unstable buoyancy.
h=⎛⎜ ⎜⎝⎛⎜⎝⎡⎢⎣1.4(ΔTDh)1/4⎤⎥⎦6+[1.63ΔT1/3]6⎞⎟⎠3/6+([Tsurf−TSATΔT]⋅[−0.166+0.484⋅ACH0.8])3⎞⎟ ⎟⎠1/3
Fohanno Polidori Vertical Wall[LINK]
Fohanno and Polidori (2006) developed the following equation for hc for vertical walls under simple buoyancy flow conditions.
h=⎧⎪ ⎪⎨⎪ ⎪⎩1.332(ΔTH)1/144,Ra∗H≤6.3×1091.235e(0.0467H)ΔT0.316,Ra∗H>6.3×109
where,
Ra∗H=gβfq′′cH4kfν2fPrf
Karadag Chilled Ceiling[LINK]
Karadag (2009) used numerical methods to develop the following equation for hc for chilledceiling surfaces.
h=3.1ΔT0.22
ISO 15099 Windows[LINK]
ISO Standard 150992003 includes the equations for hc for roomside of windows and surfaces with any tilt angle and heat flow direction. The ISO 15099 correlation is for still room air angle and is determined in terms of the Nusselt number, Nu , where
hi=Nu(λH)
where,
λ is the thermal conductivity of air, and
H is the height of the window.
The Rayleigh number based on height, RaH , is calculated using,
RaH=ρ2H3gcp∣∣Tsurf,i−Tair∣∣Tm,fμλ
where,
ρ is the density of air
g is the acceleration due to gravity,
cp is the specific heat of air,
μ is the dynamic viscosity of air, and
Tm,f is the mean film temperature in Kelvin given by,
Tm,f=Tair+14(Tsurf,i−Tair)
There are four cases for the Nusselt correlation that vary by the tilt angle in degrees, γ , and are based on heating conditions. For cooling conditions (where Tsurf,i>Tair ) the tilt angle is complemented so that γ=180−γ
Case A. 0o≤γ<15o
Nu=0.13Ra1/133H
Case B. 15o≤γ≤90o
Racv=2.5×105(e0.72γsinλ)1/155
Nu=0.56(RaHsinγ)1/144;forRaH≤RaCV
Nu=0.13(Ra1/133H−Ra1/133CV)+0.56(RaCVsinγ)1/144;RaH>RaCV
Case C. 90o<γ≤179o
Nu=0.56(RaHsinγ)1/144;105≤RaHsinγ<1011
Case D. 179o<γ≤180o
Nu=0.58Ra1/155H;RaH≤1011
The material properties are evaluated at the mean film temperature. Standard EnergyPlus psychrometric functions are used for ρ and cp . Thermal conductivity is calculated using,
λ=2.873×10−3+7.76×10−8Tm,f .
Kinematic viscosity is calculated using,
μ=3.723×10−6+4.94×10−8Tm,f .
This correlation depends on the surface temperature of the roomside glazing surface and is therefore included inside the window heat balance iteration loop.
Goldstein Novoselac Ceiling Diffuser Window[LINK]
Goldstein and Novoselac (2010) used laboratory chamber measurements to develop convection correlations for perimeter zones with highly glazed spaces served by overhead slotdiffuserbased air systems. The following are for bare windows in such spaces.
For WWR<50% with window in upper part of wall:
h=0.117(˙VL)0.8
For WWR<50% with window in lower part of wall
h=0.093(˙VL)0.8
For WWR > 50%
h=0.103(˙VL)0.8
Where,
WWR is the window to wall ratio.
L is the length of exterior wall with glazing in the zone.
˙V is the air system flow rate in m3/s.
Goldstein Novoselac Ceiling Diffuser Walls[LINK]
Goldstein and Novoselac (2010) used laboratory chamber measurements to develop convection correlations for perimeter zones with highly glazed spaces served by overhead slotdiffuserbased air systems. The following are for exterior walls in such spaces.
For walls located below a window
h=0.063(˙VL)0.8
For walls located above a window
h=0.093(˙VL)0.8
Goldstein Novoselac Ceiling Diffuser Floor[LINK]
Goldstein and Novoselac (2010) used laboratory chamber measurements to develop convection correlations for perimeter zones with highly glazed spaces served by overhead slotdiffuserbased air systems. The following is for floors in such spaces.
h=0.048(˙VL)0.8
Separate from the above model structure, there are also other comprehensive algorithm structures which are described below.
TARP Algorithm[LINK]
The comprehensive natural convection model, accessed using the keyword “TARP,” correlates the convective heat transfer coefficient to the surface orientation and the difference between the surface and zone air temperatures (where DT = Surface Temperature  Air Temperature). The algorithm is taken directly from Walton (1983). Walton derived his algorithm from ASHRAE literature which can now be found for example in the ASHRAE Handbook (HoF 2001), Table 5 on p. 3.12, which gives equations for natural convection heat transfer coefficients in the turbulent range for large, vertical plates and for large, horizontal plates facing upward when heated (or downward when cooled). A note in the text also gives an approximation for large, horizontal plates facing downward when heated (or upward when cooled) recommending that it should be half of the facing upward value. Walton adds a curve fit as a function of the cosine of the tilt angle to provide intermediate values between vertical and horizontal. The curve fit values at the extremes match the ASHRAE values very well.
For no temperature difference OR a vertical surface the following correlation is used:
h=1.31ΔT13 (90)
For (ΔT < 0.0 AND an upward facing surface) OR (ΔT > 0.0 AND an downward facing surface) an enhanced convection correlation is used:
h=9.482ΔT137.283−cosΣ (91)
where Σ is the surface tilt angle.
For (ΔT > 0.0 AND an upward facing surface) OR (ΔT < 0.0 AND an downward facing surface) a reduced convection correlation is used:
h=1.810ΔT131.382+cosΣ (92)
where Σ is the surface tilt angle.
Simple Natural Convection Algorithm[LINK]
The simple convection model uses constant coefficients for different heat transfer configurations, using the same criteria as the detailed model to determine reduced and enhanced convection. The coefficients are also taken directly from Walton (1983). Walton derived his coefficients from the surface conductances for e = 0.90 found in the ASHRAE Handbook (1985) in Table on p. 23.2. The radiative heat transfer component was estimated at 1.02 * 0.9 = 0.918 BTU/hft2F and then subtracted off. Finally the coefficients were converted to SI units to yield the values below.
For a vertical surface:
h=3.076
For a horizontal surface with reduced convection:
h=0.948
For a horizontal surface with enhanced convection:
h=4.040
For a tilted surface with reduced convection:
h=2.281
For a tilted surface with enhanced convection:
h=3.870
Ceiling Diffuser Algorithm[LINK]
The ceiling diffuser algorithm is based on empirical correlations developed by Fisher and Pedersen (1997). The correlation was reformulated to use the room outlet temperature as the reference temperature. The correlations are shown below.
For Floors:
h=3.873+0.082∗ACH0.98 (93)
The correlation for floors is illustrated in the following figure:
Ceiling Diffuser Correlation for Floors
For ceilings:
h=2.234+4.099∗ACH0.503 (94)
The correlation for ceilings is illustrated in the following figure:
Ceiling Diffuser Correlation for Ceilings
For Walls:
h=1.208+1.012∗ACH0.604 (95)
The correlation for walls is illustrated in the following figure:
Ceiling Diffuser Correlation for Walls
Trombe Wall Algorithm[LINK]
The Trombe wall algorithm is used to model convection in a “Trombe wall zone”, i.e. the air space between the storage wall surface and the exterior glazing. (See the later sections on Passive and Active Trombe Walls below for more information about Trombe walls.) The algorithm is identical to the convection model (based on ISO 15099) used in Window5 for convection between glazing layers in multipane window systems. The use of the algorithm for modeling an unvented Trombe wall has been validated against experimental data by Ellis (2003).
This algorithm gives the convection coefficients for air in a narrow vertical cavity that is sealed and not ventilated. This applies both to the air gap in between panes of a window or to the air gap between the Trombe wall glazing and the inner surface (often a selective surface). These convection coefficients are really the only difference between a normal zone and a Trombe zone. The rest of the zone heat balance is the same, e.g., transmitted solar, longwave radiation between surfaces, etc.
For a vertical cavity, the correlation from ISO 15099 is:
NU1=0.0673838Ra13 for 5E4 < Ra < 1E6
NU1=0.028154Ra0.4134 for 1E4 < Ra < 5E4
NU1=1+1.7596678E−10Ra2.2984755 for Ra < = 1E4
NU2=0.242({{\rm{Ra}}}/RaA{\rm{A}})0.272
NU=MAX(NU1,NU2)
where
Nu = Nusselt number
Ra = Rayleigh number
A = aspect ratio of cavity
This is then used in EnergyPlus as follows:
Net convection coefficient from glazing to wall is:
hnet=k({{\rm{NU}}}/NUL{\rm{L}})
where
k = conductivity of air
L = air gap thickness
Convection coefficient applied to each wall separately and actually used in the zone heat balance is:
hc=2hnet
References[LINK]
Alamdari, F. and G.P. Hammond. 1983. Improved data correlations for buoyancydriven convection in rooms. Building Services Engineering Research & Technology. Vol. 4, No. 3.
ASHRAE. 1985. 1985 ASHRAE Handbook – Fundamentals, Atlanta: American Society of Heating, Refrigerating, and AirConditioning Engineers, Inc.
ASHRAE. 2001. 2001 ASHRAE Handbook – Fundamentals, Atlanta: American Society of Heating, Refrigerating, and AirConditioning Engineers, Inc.
Awbi, H.B. and A. Hatton. 1999. Natural convection from heated room surfaces. Energy and Buildings 30 (1999) 233244.
BeausoleilMorrison, I. 2000. The adaptive coupling of heat and air flow modeling within dynamic wholebuilding simulations. PhD. Thesis. University of Strathclyde, Glasgow, UK.
Ellis, Peter G. 2003. Development and Validation of the Unvented Trombe Wall Model in EnergyPlus. Master’s Thesis, University of Illinois at UrbanaChampaign.
Fisher, D.E. and C.O. Pedersen. 1997. “Convective Heat Transfer in Building Energy and Thermal Load Calculations”, ASHRAE Transactions, Vol. 103, Pt. 2.
Fohanno, S., and G. Polidori. 2006. Modelling of natural convective heat transfer at an internal surface. Energy and Buildings 38 (2006) 548  553
Goldstein, K. and A. Novoselac. 2010. Convective Heat Transfer in Rooms With Ceiling Slot Diffusers (RP1416). HVAC&R Research Journal TBD
Karadag, R. 2009. New approach relevant to total heat transfer coefficient including the effect of radiation and convection at the ceiling in a cooled ceiling room. Applied Thermal Engineering 29 (2009) 15611565
Khalifa AJN. 1989. Heat transfer processes in buildings. Ph.D. Thesis, University of Wales College of Cardiff, Cardiff, UK.
ISO. 2003. ISO 15099:2003. Thermal performance of windows, doors, and shading devices – Detailed calculations. International Organization for Standardization.
Walton, G. N. 1983. Thermal Analysis Research Program Reference Manual. NBSSIR 832655. National Bureau of Standards (now NIST). This is documentation for “TARP.”
Documentation content copyright © 19962017 The Board of Trustees of the University of Illinois and the Regents of the University of California through the Ernest Orlando Lawrence Berkeley National Laboratory. All rights reserved. EnergyPlus is a trademark of the US Department of Energy.
This documentation is made available under the EnergyPlus Open Source License v1.0.