Solar Collectors[LINK]

Solar collectors are devices that convert solar energy into thermal energy by raising the temperature of a circulating heat transfer fluid. The fluid can then be used to heat water for domestic hot water usage or space heating. Flat-plate solar collectors using water as the heat transfer fluid, Integral-Collector Storage solar collectors using water and unglazed transpired solar collectors using air are currently the only types of collector available in EnergyPlus.

Flat-Plate Solar Collectors[LINK]

The input object SolarCollector:FlatPlate:Water provides a model for flat-plate solar collectors that are the most common type of collector. Standards have been established by ASHRAE for the performance testing of these collectors (ASHRAE 1989; 1991) and the Solar Rating and Certification Corporation (SRCC) publishes a directory of commercially available collectors in North America (SRCC 2003).

The EnergyPlus model is based on the equations found in the ASHRAE standards and Duffie and Beckman (1991). This model applies to glazed and unglazed flat-plate collectors, as well as banks of tubular, i.e. evacuated tube, collectors.

Solar and Shading Calculations[LINK]

The solar collector object uses a standard EnergyPlus surface in order to take advantage of the detailed solar and shading calculations. Solar radiation incident on the surface includes beam and diffuse radiation, as well as radiation reflected from the ground and adjacent surfaces. Shading of the collector by other surfaces, such as nearby buildings or trees, is also taken into account. Likewise, the collector surface can shade other surfaces, for example, reducing the incident radiation on the roof beneath it.

Thermal Performance[LINK]

The thermal efficiency of a collector is defined as the ratio of the useful heat gain of the collector fluid versus the total incident solar radiation on the gross surface area of the collector.

$$\eta =\frac{\left(q/A\right)}{I_{solar}}$$

where:

q is the useful heat gain

A is the gross area of the collector

$I_{solar}$ is the total incident solar radiation.

Notice that the efficiency h is only defined for I${}_{solar}$ > 0.

An energy balance on a solar collector with double glazing shows relationships between the glazing properties, absorber plate properties, and environmental conditions.

$$\frac{q}{A}=I_{solar}{\tau }_{g1}{\tau }_{g2}{\alpha }_{abs}-\frac{T_{abs}^4-T_{g2}^4}{R_{rad}}-\frac{T_{abs}-T_{g2}}{R_{conv}}-\frac{T_{abs}-T_{air}}{R_{cond}}$$

where:

${\tau }_{g1}$ is the transmittance of the first glazing layer

${\tau }_{g2}$ is the transmittance of the second glazing layer

${\alpha }_{abs}$ is the absorptance of the absorber plate

$R_{rad}$ is the radiative resistance from absorber to inside glazing

$R_{conv}$ is the convective resistance from absorber to inside glazing

$R_{cond}$ is the conductive resistance from absorber to outdoor air through the insulation

$T_{abs}$ is the temperature of the absorber plate

$T_{g2}$ is the temperature of the inside glazing

$T_{air}$ is the temperature of the outdoor air.

The equation above can be approximated with a simpler formulation as:

$$\frac{q}{A}=F_R\left[I_{solar}\left(\tau \alpha \right)-U_L(T_{in}-T_{air})\right]$$

where:

$F_R$ is an empirically determined correction factor

($\tau \alpha$) is the product of all transmittance and absorptance terms

$U_L$ is the overall heat loss coefficient combining radiation, convection, and conduction terms

$T_{in}$ is the inlet temperature of the working fluid.

Substituting this into Equation [eq:SolarCollectorEta],

$$\eta =F_R\left(\tau \alpha \right)-F_R{U}_L\frac{(T_{in}-T_{air})}{I_{solar}}$$

A linear correlation can be constructed by treating F${}_R$($\tau \alpha$) and -F${}_R$U${}_L$ as characteristic constants of the solar collector:

$$\eta =c_0+c_1\frac{\left(T_{in}-T_{air}\right)}{I_{solar}}$$

Similarly, a quadratic correlation can be constructed using the form:

$$\eta =c_0+c_1\frac{\left(T_{in}-T_{air}\right)}{I_{solar}}+c_2\frac{{\left(T_{in}-T_{air}\right)}^2}{I_{solar}}$$

Both first- and second-order efficiency equation coefficients are listed in the Directory of SRCC Certified Solar Collector Ratings.

Incident Angle Modifiers[LINK]

As with regular windows the transmittance of the collector glazing varies with the incidence angle of radiation. Usually the transmittance is highest when the incident radiation is normal to the glazing surface. Test conditions determine the efficiency coefficients for normal incidence. For off-normal angles, the transmittance of the glazing is modified by an incident angle modifier coefficient.

$$K_{\tau \alpha }=\frac{\left(\tau \alpha \right)}{{\left(\tau \alpha \right)}_n}$$

Additional testing determines the incident angle modifier as a function of incident angle $\theta$. This relationship can be fit to a first-order, linear correlation:

$$K_{\tau \alpha }=1+b_0\left(\frac{1}{\cos \theta }-1\right)$$

or a second-order, quadratic correlation:

$$K_{\tau \alpha }=1+b_0\left(\frac{1}{\cos \theta }-1\right)+b_1{\left(\frac{1}{\cos \theta }-1\right)}^2$$

The incident angle modifier coefficients b${}_0$ and b${}_1$ are usually negative, although some collectors have a positive value for b${}_0$. Both first- and second-order incident angle modifier equation coefficients are listed in the Directory of SRCC Certified Solar Collector Ratings.

The SRCC incident angle modifier equation coefficients are only valid for incident angles of 60 degrees or less. Because these curves can be valid yet behave poorly for angles greater than 60 degree, the EnergyPlus model cuts off collector gains for incident angles greater than 60 degrees.

For flat-plate collectors, the incident angle modifier is generally symmetrical. However, for tubular collectors the incident angle modifier is different depending on whether the incident angle is parallel or perpendicular to the tubes. These are called bi-axial modifiers. Some special flat-plate collectors may also exhibit this asymmetry. The current model cannot yet handle two sets of incident angle modifiers. In the meantime it is recommended that tubular collectors be approximated with caution using either the parallel or perpendicular correlation.

Incident angle modifiers are calculated separately for sun, sky, and ground radiation. The net incident angle modifier for all incident radiation is calculated by weighting each component by the corresponding modifier.

$$K_{\tau \alpha ,net}=\frac{I_{beam}{K}_{\tau \alpha ,beam}+I_{sky}{K}_{\tau \alpha ,sky}+I_{gnd}{K}_{\tau \alpha ,gnd}}{I_{beam}+I_{sky}+I_{gnd}}$$

For sky and ground radiation the incident angle is approximated using Brandemuehl and Beckman’s equations:

$${\theta }_{sky}=59.68-0.1388\varphi +0.001497{\varphi }^2$$

$${\theta }_{ground}=90.0-0.5788\varphi +0.002693{\varphi }^2$$

where $\varphi$ is the surface tilt in degrees.

The net incident angle modifier is then inserted into the useful heat gain Equation [eq:SolarCollectorqOverA697]:

$$\frac{q}{A}=F_R\left[I_{solar}{K}_{\tau \alpha ,net}{\left(\tau \alpha \right)}_n-U_L(T_{in}-T_{air})\right]$$

Equation [eq:SolarCollectorEta698] is also modified accordingly.

$$\eta =F_R{K}_{\tau \alpha ,net}(\tau \alpha {)}_n-F_R{U}_L\frac{(T_{in}-T_{air})}{I_{solar}}$$

Outlet Temperature[LINK]

Outlet temperature is calculated using the useful heat gain q as determined by Equation [eq:SolarCollectorqOverA707], the inlet fluid temperature T${}_{in}$, and the mass flow rate available from the plant simulation:

$$\frac{q}{A}=\dot{m}{c}_p(T_{out}-T_{in})$$

where:

$\dot{m}$ is the fluid mass flow rate through the collector

$c_p$ is the specific heat of the working fluid.

Solving for T${}_{out}$,

$$T_{out}=T_{in}+\frac{q}{\dot{m}{c}_{p}A}$$

If there is no flow through the collector, T${}_{out}$ is the stagnation temperature of the fluid. This is calculated by setting the left side of Equation [eq:SolarCollectorqOverA707] to zero and solving for T${}_{in}$ (which also equals T${}_{out}$ for the no flow case).

References[LINK]

ASHRAE. 1989. ASHRAE Standard 96-1980 (RA 89): Methods of Testing to Determine the Thermal Performance of Unglazed Flat-Plate Liquid-Type Solar Collectors. Atlanta: American Society of Heating, Refrigerating and Air-Conditioning Engineers, Inc.

ASHRAE. 1991. ASHRAE Standard 93-1986 (RA 91): Methods of Testing to Determine the Thermal Performance of Solar Collectors. Atlanta: American Society of Heating, Refrigerating and Air-Conditioning Engineers, Inc.

Duffie, J. A., and Beckman, W. A. 1991. Solar Engineering of Thermal Processes, Second Edition. New York: Wiley-Interscience.

Solar Rating and Certification Corporation. 2004. Directory of SRCC Certified Solar Collector Ratings, OG 100. Cocoa, Florida: Solar Rating and Certification Corporation.

Integral-collector-storage (ICS) Solar Collector[LINK]

Solar collectors with integral storage unit models use SolarCollector:IntegralCollectorStorage object, and the characteristics parameter inputs of this collector are provided by the SolarCollectorPerformance:IntegralCollectorStorage object. This model is based on detailed Energy Balance equations of solar collectors that integrates storage in it. This model has two options to represent the collector bottom outside boundary conditions: AmbientAir, and OtherSideConditionsModel. AmbientAir simply applies outside air temperature using combined convection and radiation conductance, and the OtherSideConditionsModel applies combined radiation and convection models that exist in a naturally ventilated cavity to represent the collector bottom outside boundary condition. The later boundary condition accounts for the shading of the collector on the underlying surface, hence, the ICS collector can be assumed as an integral part of the building envelope. Schematic diagram of a rectangular ICS solar collector is shown in Figure 1 below:

Solar and Shading Calculations[LINK]

The solar collector object uses a standard EnergyPlus surface in order to take advantage of the detailed solar and shading calculations. Solar radiation incident on the surface includes beam and diffuse radiation, as well as radiation reflected from the ground and adjacent surfaces. Shading of the collector by other surfaces, such as nearby buildings or trees, is also taken into account. Likewise, the collector surface shades the roof surface beneath it, hence no direct solar radiation incident on the roof surface. The collector and the roof outside boundary conditions should be specified as OtherSideConditionModel to account for solar collector shading impact on the roof surface.

Mathematical Model[LINK]

The integral-collector-storage (ICS) solar collector is represented using two transient energy balance equations shown below. These equations represent the energy balance equation for the absorber plate, and the water in the collector.

$$m_p{C}_p\frac{d{T}_p}{dt}=A\cdot \left[{\left(\tau \alpha \right)}_e\cdot I_t-h_{pw}\left(T_p-T_w\right)-U_t\left(T_p-T_a\right)\right]$$

$$m_w{C}_w\frac{d{T}_w}{dt}=A\cdot \left[h_{pw}\left(T_p-T_w\right)-U_b\left(T_w-T_{OSC}\right)-U_s\left(T_w-T_a\right)\right]-{\dot{m}}_w{C}_w\left(T_w-T_{wi}\right)$$

where:

$m_p{C}_p$ is the thermal capacity of the absorber surface, J/${}^{\circ }$C

$A$ is the collector gross area, m${}^2$

${\left(\tau \alpha \right)}_e$ is the transmittance-absorptance product of the absorber plate and cover system

$I_t$ is the total solar irradiation (W/m${}^2$)

$h_{pw}$ is the convective heat transfer coefficient from absorber plate to water (W/m${}^2$-K)

$m_a$ is the overall heat loss coefficient from absorber to the ambient air (W/m${}^2$-K)

$T_p$ is the absorber plate average temperature (${}^{\circ }$C)

$T_w$ is the collector water average temperature (${}^{\circ }$C)

$T_a$ is the ambient air temperature (${}^{\circ }$C)

$m_w{C}_{pw}$ is the thermal capacity of the water mass in the collector (J/${}^{\circ }$C)

$U_s$ is the area-weighted conductance of the collector side insulation (W/m${}^2$-K)

$U_b$ is the conductance of the collector bottom insulation (W/m${}^2$-K)

$T_{osc}$ is the outside temperature of bottom insulation determined from the other side condition model (${}^{\circ }$C)

$T_{wi}$ is the entering makeup or mains water temperature (${}^{\circ }$C)

${\dot{m}}_w{C}_w$ is the water capacity flow through the collector (W/${}^{\circ }$C).

The other side condition model boundary condition represented by the T${}_{osc}$, allows us to apply a realistic outside boundary condition for a collector mounted on a building roof. This also accounts for the shading impact of the collector on the under-laying surface (roof). On the other hand if ambient air boundary condition is specified, then the collector does not shade the underlying surface it is mounted on.

The two energy balance equation can be written as non-homogeneous first order DE with constant coefficients. The initial conditions for these equations are the absorber plate average temperature and the collector water average temperature at previous time steps.

$$\frac{d{T}_p}{dt}=a_1{T}_p+a_2{T}_w+a_3$$

$$\frac{d{T}_w}{dt}=b_1{T}_p+b_2{T}_w+b_3$$

$$a_1=-\left(A\cdot h_{pw}+A\cdot U_t\right)/\left(m_p\cdot C_p\right)$$

$$a_2=A\cdot h_{pw}{T}_w/\left(m_p\cdot C_p\right)$$

$$a_3=A\cdot {\left(\tau \alpha \right)}_e\cdot I_t+A\cdot U_t{T}_a$$

$$b_1=A\cdot h_{pw}{T}_p/\left(m_w\cdot C_w\right)$$

$$b_2=-\left(A\cdot h_{pw}+A\cdot U_b+A\cdot U_s+{\dot{m}}_w{C}_w\right)$$

$$b_3=\left(A\cdot U_b{T}_{osc}+A\cdot U_s{T}_a+{\dot{m}}_w{C}_w{T}_{wi}\right)$$

The two coupled first order differential equation are solved analytically. Auxiliary equation of the the coupled homogeneous differential equation is given by:

$${\lambda }^2-(a_1+b_2)+(a_1{b}_2-a_2{b}_1)=0$$

This auxiliary quadratic equation has always two distinct real roots (${\lambda }_1$ and ${\lambda }_2$) hence the solution of the homogeneous equation is exponential, and the general solutions of the differential equations are given by:

$$T_p=c_1{e}^{{\lambda }_{1}t}+c_2{e}^{{\lambda }_{2}t}+A$$

$$T_w=r_1{c}_1{e}^{{\lambda }_{1}t}+r_2{c}_2{e}^{{\lambda }_{2}t}+B$$

The constant terms A and B are the particular solution of the non-homogeneous differential equations, the coefficients of the exponential terms (c${}_1$, c${}_2$, r${}_1$, and r${}_2$) are determined from the initial conditions of the absorber and collector water temperatures (T${}_{p0}$, T${}_{w0}$) and are given by:

$$r_1=({\lambda }_1-a_1)/a_2;r_2=({\lambda }_2-a_1)/a_2$$

$$A=(-a_3{b}_2+b_3{a}_2)/(a_1{b}_2-b_1{a}_2);B=(-a_1{b}_3+b_1{a}_3)/(a_1{b}_2-b_1{a}_2)$$

$$c_1=(r_2{T}_{p0}-T_{w0}-r_{2}A+B)/(r_2-r_1);c_2=(T_{w0}-r_1{T}_{p0}+r_{1}A-B)/(r_2-r_1)$$

Thermal Network Model:[LINK]

The thermal network model requires energy balance for each of the collector covers as well. The heat balance equation of the collector covers is assumed to obey steady state formulation by ignoring their thermal mass. The thermal-network representation of the ICS collector is shown in Figure 2. Also, the heat balance at each cover surface requires knowledge of the amount of solar fraction absorbed, which is determined from the ray tracing analysis. For the thermal network model shown above the overall top heat loss coefficient is determined from combination of the resistances in series as follows:

$$U_t={\left[R_1+R_2+R_3\right]}^{-1}$$

or

$$U_t={\left[\frac{1}{h_{c,c1-a}+h_{r,c1-a}}+\frac{1}{h_{c,c2-c1}+h_{r,c2-c1}}+\frac{1}{h_{c,p-c2}+h_{r,p-c2}}\right]}^{-1}$$

The convection and radiation heat transfer coefficients in equation above are calculated based on temperatures at the previous time step and determined as described in the Heat Transfer Coefficients section.

Collector Cover Heat Balance

Ignoring the thermal mass of the collector cover, a steady state heat balance equations are formulated for each cover that allow us to determine cover temperatures. The cover surface heat balance representation is shown in Figure 3 below.

The steady state cover heat balance equation is given by:

$$q_{LWR,1}+q_{CONV,1}+q_{solar,abs}+q_{LWR,2}+q_{CONV,2}=0$$

Linearizing the longwave radiation exchange and representing the convection terms using the classical equation for Newton’s law of cooling, the equations for the temperatures of covers 1 and 2 are given by:

$$T_{c1}=\frac{{\alpha }_{c1}{I}_t+h_{r,c1-a}{T}_a+h_{c,c1-a}{T}_a+h_{r,c2-c1}{T}_{c2}+h_{c,c2-c1}{T}_{c2}}{h_{r,c1-a}+h_{c,c1-a}+h_{r,c2-c1}+h_{c,c2-c1}}$$

$$T_{c2}=\frac{{\alpha }_{c2}{I}_t+h_{r,c2-c1}{T}_{c1}+h_{c,c2-c1}{T}_{c1}+h_{r,p-c2}{T}_p+h_{c,p-c2}{T}_p}{h_{r,c2-c1}+h_{c,c2-c1}+h_{r,p-c2}+h_{c,p-c2}}$$

where:

${\alpha }_c$ is the weighted average solar absorptance of covers 1 and 2

$h_{r,c1-a}$ is the adjusted radiation heat transfer coefficient between cover 1 and the ambient air (W/m${}^2$-K)

$h_{c,c1-a}$ is the convection heat transfer coefficient between cover 1 and the ambient (W/m${}^2$-K)

$h_{r,c2-c1}$ is the radiation heat transfer coefficient between covers 1 and 2 (W/m${}^2$-K)

$h_{c,c2-c1}$ is the convection heat transfer coefficient between covers 1 and 2 (W/m${}^2$-K)

$h_{r,p-c2}$ is the radiation heat transfer coefficient between covers 2 and the absorber plate (W/m${}^2$-K)

$h_{c,p-c2}$ is the convection heat transfer coefficient between covers 2 and the absorber plate (W/m${}^2$-K)

$q_{LWR,1}$ is the longwave radiation exchange flux on side 1 of the collector cover (W/m${}^2$)

$q_{CONV,1}$ is the convection heat flux on side 1 of the collector cover (W/m${}^2$)

$q_{LWR,2}$ is the longwave radiation exchange flux on side 2 of the collector cover (W/m${}^2$)

$q_{CONV,2}$ is the convection heat flux on side 2 of the collector cover (W/m${}^2$)

$q_{solar,abs}$ is the net solar radiation absorbed by the collector cover (W/m${}^2$)

$R$ is the thermal resistance for each section along the heat flow path (m${}^2$K/W).

Other Side Condition Model[LINK]

ICS Solar Collectors are commonly mounted on building heat transfer surfaces hence the collectors shade the underlying heat transfer surface and require a unique boundary condition that reflects the air cavity environment created between the bottom of the collector surface and the underlying surface. The other side condition model that allows us to estimate the other side temperature, T${}_{osc}$, can be determined based on steady state heat balance using the known collector water temperature at the previous time step.

Ignoring thermal mass of the collector bottom insulation, steady state surface heat balance can be formulated on the outer plane of the collector bottom surface facing the cavity as shown in Figure 4. The heat balance equation on the outer plane of the collector bottom surface is given by:

$$q_{cond}+q_{conv,cav}+q_{rad,cav}=0$$

Substituting the equations for each term in the above equation yields:

$$U_L(T_w-T_{osc})+h_{c,cav}(T_{a,cav}-T_{osc})+h_{r,cav}(T_{so}-T_{osc})=0$$

Simplifying yields the bottom insulation other side condition temperature:

$$T_{osc}=\frac{U_L{T}_w+h_{c,cav}{T}_{a,cav}+h_{r,cav}{T}_{so}}{U_L+h_{c,cav}+h_{r,cav}}$$

The cavity air temperature is determined from cavity air heat balance as follows:

$$T_{a,cav}=\frac{h_{c,cav}A{T}_{osc}+{\dot{m}}_{vent}Cp{T}_a+h_{c,cav}A{T}_{so}}{h_{c,cav}A+{\dot{m}}_{vent}Cp+h_{c,cav}A}$$

where:

$h_{r,cav}$ is the linearized radiation coefficient for underlying surface in the cavity (W/m${}^2$-K)

$h_{c,cav}$ is the convection coefficient for underlying surface in the cavity (W/m${}^2$-K)

$T_{so}$ is the outside face temperature of the underlying heat transfer surface (${}^{\circ }$C)

${\dot{m}}_{vent}$ is the air mass flow rate due to natural ventilation (kg/s)

$q_{cond}$ is the conduction heat flux though the insulation and bottom (W/m${}^2$)

$q_{conv,cav}$ is the convection heat flux between the collector bottom outside surface and the cavity air (W/m${}^2$)

$q_{rad,cav}$ is the longwave radiation exchange flux between the collector bottom outside surface and the outside surface of the underlying surface (W/m${}^2$).

The cavity air temperature is determined from the cavity air energy balance. The air heat balance requires the ventilated cavity air natural ventilation rates. The calculation of the ventilation rate is described else where in this document. The SurfaceProperty:ExteriorNaturalVentedCavity, object is required to describe the surface properties, the characteristics of the cavity and opening for natural ventilation.

Heat Transfer Coefficients[LINK]

The equations used to determine for the various heat transfer coefficients in the absorber and water heat balance equations are given below. The absorbed solar energy is transferred to the water by convection. Assuming natural convection dominated heat transfer for a hot surface facing down and a clod surface facing down the following correlation for Nusselt number by Fujii and Imura (1972). The Nusselt number for hot surface facing down ward is given by:

$$Nu=0.56{\left(Gr\cdot Pr\cdot cos\theta \right)}^{1/5}{10}^5<Gr\cdot Pr<{10}^{11}$$

The Nusselt number for hot surface upward and cold surface facing down is given by:

$$Nu=0.13{\left(Gr\cdot Pr\right)}^{1/3}Gr\cdot Pr<5.0\times {10}^8$$

$$Nu=0.16{\left(Gr\cdot Pr\right)}^{1/3}Gr\cdot Pr>5.0\times {10}^8$$

$$Gr=g{\beta }_v\left(T_p-T_w\right)L_c^3/g{\beta }_v\left(T_p-T_w\right)L_c^3{\nu }^2{\nu }^2$$

$$Pr=\nu /\nu \alpha \alpha$$

$$T_r=T_p-0.25\left(T_p-T_w\right)$$

$$h_w=Nu\cdot k/Nu\cdot k{L}_c{L}_c$$

where:

$\theta$ is the angle of inclination of the collector to the vertical (radians)

g is the gravitation force constant, 9.806 (m/s${}^2$)

$Q_n$ is the reference properties where the thermo-physical properties are calculated (${}^{\circ }$C)

$L_c$ is the characteristic length for the absorber plate (m)

k is the thermal conductivity of water at reference temperature (W/m-K)

$\nu$ is the kinematic viscosity of water at reference temperature (m${}^2$/s)

$\alpha$ is the thermal diffusivity of water at reference temperature (m${}^2$/s)

${\beta }_v$ is the volumetric expansion coefficient (1/K) evaluated at Tv, Tv = Tw+0.25(Tp-Tw)

Nu is the Nusselt number calculated for water properties at the reference temperature

Gr is the Grashof number calculated for water properties at the reference temperature

Pr is the Prandtl number calculated for water properties at the reference temperature.

The various radiation and convection heat transfer coefficients are given by the following equations. The convection heat transfer coefficients between the covers and the absorber plate are estimated from the empirical correlation for the Nusselt number for air gap between two parallel plates developed by Hollands et al. (1976) is:

$$N{u}_a=1+1.44\left\{1-\frac{1708{\left(\sin 1.8\beta \right)}^{1.6}}{Ra\cos \beta }\right\}\left\{1-\frac{1708}{Ra\cos \beta }\right\}+\left\{{\left(\frac{Ra\cos \beta }{5830}\right)}^{\text{1}/13\text{3}}-1\right\}$$

$$h_c=Nu\cdot k/Nu\cdot kLL$$

$$h_{rp-c2}=\frac{\sigma \left(T_p+T_{c2}\right)\left(T_p^2+T_{c2}^2\right)}{1/{\epsilon }_p+1/{\epsilon }_{c2}-1}$$

$$h_{rc1-c2}=\frac{\sigma \left(T_{c1}+T_{c2}\right)\left(T_{c1}^2+T_{c2}^2\right)}{1/{\epsilon }_{c1}+1/{\epsilon }_{c2}-1}$$

The long wave radiation exchange coefficients between the outer collector cover and the sky and ground referencing the ambient air temperature for mathematical simplification are given.

$$h_{rc1-s}=\frac{F_s{\epsilon }_{c1}\sigma \left(T_{c1}+T_s\right)\left(T_{c1}^2+T_s^2\right)\left(T_{c1}-T_s\right)}{\left(T_{c1}-T_a\right)}$$

$$h_{rc1-g}=\frac{F_g{\epsilon }_{c1}\sigma \left(T_{c1}+T_g\right)\left(T_{c1}^2+T_g^2\right)\left(T_{c1}-T_g\right)}{\left(T_{c1}-T_a\right)}$$

$$h_{cc1-a}=h_{cc1-s}+h_{cc1-g}$$

The convection heat transfer coefficient from the outer cover to the surrounding air is given by:

$$h_{cc1-a}=2.8+3.0V_w$$

When the bottom surface boundary condition is AmbientAir, the combined conductance from the outer cover to the surrounding is calculated from the equation below (Duffie and Beckman, 1991).

$$h_{comb}=5.7+3.8V_w$$

The overall loss coefficient through the bottom and side of the collector-storage is estimated as follows:

$$U_b=U_{Lb}\left(A_b/A\right)$$

$$U_s={\left[\frac{1}{U_{Ls}\left(A_s/A\right)}+\frac{1}{h_{comb}}\right]}^{-1}$$

where:

${\epsilon }_{c1}$ is the thermal emissivity of collector cover 1

${\epsilon }_{c2}$ is the thermal emissivity of collector cover 2

$F_s$ is the view factor from the collector to the sky

$F_g$ is the view factor from the collector to the ground

$T_{c1}$ is the temperature of collector cover 1, (K)

$T_{c2}$ is the temperature of collector cover 2, (K)

$T_s$ is the sky temperature, (K)

$T_g$ is the ground temperature, (K)

k is the thermal conductivity of air (W/m K)

L is the air gap between the covers (m)

$\beta$ is the inclination of the plates or covers to the horizontal (radians)

$V_w$ is the wind speed (m/s)

$U_{Lb}$ is the user specified bottom heat loss conductance (W/m${}^2$-K)

$U_{Ls}$ is the user specified side heat loss conductance (W/m${}^2$-K)

$A_b$ is the collector bottom heat transfer area (m${}^2$)

$A_s$ is the collector side area (m${}^2$)

$h_{comb}$ is the combined conductance from the outer cover to the ambient air (W/m${}^2$-K).

Transmittance-Absorptance Product

The transmittance-absorptance product of solar collector is determined using ray tracing method for any incident angle (Duffie and Beckman, 1991). This requires optical properties of the cover and absorber materials and the the transmittance-absorptance product for any incident angle is given by:

$${\left(\tau \alpha \right)}_{\theta }=\frac{\tau \alpha }{1-\left(1-\alpha \right){\rho }_d}$$

The transmittance of the cover system for single and two cover are given by:

$$\tau =\frac{1}{2}\left[{\left(\frac{{\tau }_1\cdot {\tau }_2}{1-{\rho }_1{\rho }_2}\right)}_{\perp }+{\left(\frac{{\tau }_1\cdot {\tau }_2}{1-{\rho }_1{\rho }_2}\right)}_{\parallel }\right]$$

$$\rho =\frac{1}{2}\left[{\left({\rho }_1+\frac{\tau \cdot {\rho }_2\cdot {\tau }_1}{{\tau }_2}\right)}_{\perp }+{\left({\rho }_1+\frac{\tau \cdot {\rho }_2\cdot {\tau }_1}{{\tau }_2}\right)}_{\parallel }\right]$$

The effective transmittance, reflectance and absorptance of a single cover are given by:

$$\tau =\frac{{\tau }_a}{2}\left\{\frac{1-r_{\perp }}{1+r_{\perp }}\left[\frac{1-r_{\perp }^2}{1-{\left(r_{\perp }{\tau }_a\right)}^2}\right]+\frac{1-r_{\parallel }}{1+r_{\parallel }}\left[\frac{1-r_{\parallel }^2}{1-{\left(r_{\parallel }{\tau }_a\right)}^2}\right]\right\}$$

$$\rho =\frac{1}{2}\left\{\left[r_{\perp }+\frac{{\left(1-r_{\perp }\right)}^2{\tau }_a^2{r}_{\perp }}{1-{\left(r_{\perp }{\tau }_a\right)}^2}\right]+\left[r_{\parallel }+\frac{{\left(1-r_{\parallel }\right)}^2{\tau }_a^2{r}_{\parallel }}{1-{\left(r_{\parallel }{\tau }_a\right)}^2}\right]\right\}$$

$$\alpha =\frac{\left(1-{\tau }_a\right)}{2}\left\{\left(\frac{1-r_{\perp }}{1-r_{\perp }{\tau }_a}\right)+\left(\frac{1-r_{\parallel }}{1-r_{\parallel }{\tau }_a}\right)\right\}$$

The transmittance of the cover system with absorption only considered t${}_a$, is defined as:

$${\tau }_a=\exp \left(-\frac{KL}{\cos {\theta }_2}\right)$$

$${\theta }_2={\sin }^{-1}\left(\sin {\theta }_1\cdot \frac{n_1}{n_2}\right)$$

The reflectance of un-polarized radiation on passing from medium 1 with reflective index n${}_1$ to medium 2 with reflective index n${}_2$ is given by:

$$r_{\perp }=\frac{{\sin }^2\left({\theta }_2-{\theta }_1\right)}{{\sin }^2\left({\theta }_2+{\theta }_1\right)}$$

$$r_{\parallel }=\frac{{\tan }^2\left({\theta }_2-{\theta }_1\right)}{{\tan }^2\left({\theta }_2+{\theta }_1\right)}$$

The sky and ground reflected diffuse radiations average equivalent incident angles are approximated by Brandemuehl and Beckman correlation (Duffie and Beckman, 1991) as follows:

$${\theta }_{sd}=59.68-0.1388\beta +0.001497{\beta }^2$$

$${\theta }_{gd}=90-0.5788\beta +0.002693{\beta }^2$$

where:

$\tau$ is the transmittance of the cover system

${\tau }_1$ is the transmittance of the cover 1

${\tau }_2$ is the transmittance of the cover 2

$\delta$ is the absorptance of the absorber plate

${\rho }_d$ is the diffuse reflectance of the inner cover

L is the thickness of a cover material (m)

K is the extinction coefficient of a cover material (m${}^{-1}$)

${\theta }_1$ is the angle of incidence (degrees)

${\theta }_2$ is the angle of refraction (degrees)

$r_{\perp }$ is the parallel component of reflected un-polarized radiation

$r_{\parallel }$ is the perpendicular component of reflected un-polarized radiation

$S_m$ is the slope of the collector (degrees)

${\theta }_{sd}$ is the equivalent incident angle for sky diffuse solar radiation (degrees)

${\theta }_{gd}$ is the equivalent incident angle for ground diffuse solar radiation (degrees).

The integral collector storage unit thermal performance parameters are calculated as follows:

$$Q_{delivered}={\dot{m}}_w{C}_w\left(T_w-T_{wi}\right)$$

$$Q_{Stored}=m_w{C}_w\frac{d{T}_w}{dt}$$

$$Q_{SkinLoss}=A\cdot \left[U_t\left(T_p-T_w\right)+U_b\left(T_w-T_{OSC}\right)+U_s\left(T_w-T_a\right)\right]$$

$${\eta }_{thermal}=\frac{m_w{C}_w\frac{d{T}_w}{dt}+{\dot{m}}_w{C}_w\left(T_w-T_{wi}\right)}{A\cdot I_t}$$

References:[LINK]

Duffie, J.A., and W.A. Beckman. 1991. Solar Engineering of Thermal Processes, 2d ed. New York: John Wiley & Sons.

Kumar, R. and M.A. Rosen. Thermal performance of integrated collector storage solar water heater with corrugated absorber surface. Applied Thermal Engineering: 30 (2010) 1764–1768.

Fujii, T., and H. Imura. Natural convection heat transfer from a plate with arbitrary inclination. International Journal of Heat and Mass Transfer: 15(4), (1972), 755-764.

Photovoltaic Thermal Flat-Plate Solar Collectors[LINK]

Photovoltaic-Thermal solar collectors (PVT) combine solar electric cells and thermal working fluid to collect both electricity and heat. Although there are currently comparatively few commercial products, PVT research has been conducted for the past 30 years and many different types of collectors have been studied. Zondag (2008) and Charalambous et. al (2007) provide reviews of the PVT literature. Because PVT is much less commercially-mature, there are no standards or rating systems such as for thermal-only, hot-water collectors. EnergyPlus currently has one simple model based on user-defined efficiencies.

The PVT models reuse the PV models for electrical production. These are described elsewhere in this document in the section on Photovoltaic Arrays-Simple Model

Simple PVT Thermal Model[LINK]

The input object SolarCollector:FlatPlate:PhotovoltaicThermal provides a simple PVT model that is provided for quick use during design or policy studies. The user simply provides values for a thermal efficiency and the incident solar heats the working fluid. The model also includes a cooling mode for air-based systems where a user-provided surface emittance is used to model cooling of the working fluid to the night sky (water-based cooling will be made available once a chilled water storage tank is available). No other details of the PVT collector’s construction are required as input data.

The simple model can heat either air or liquid. If it heats air, then the PVT is part of HVAC air system loop with air nodes connected to an air system. If it heats liquid, then the PVT is part of plant loop with nodes connected to a plant loop and the plant operating scheme determines flows.

Air-system-based PVT modeling include a modulating bypass damper arrangement. Control logic decides if the air should bypass the collector to better meet setpoint. The model requires a drybulb temperature setpoint be placed on the outlet node. The model assume the collector is intended and available for heating when the incident solar is greater than 0.3 W/m${}^2$ and otherwise it is intended for cooling. The inlet temperature is compare to the setpoint on the outlet node to determine if cooling or heating are beneficial. If they are, then the PVT thermal models are applied to condition the air stream. If they are not beneficial, then the PVT is completely bypassed and the inlet node is passed directly to the outlet node to model a completely bypassed damper arrangement. A report variable is available for bypass damper status.

Plant-based PVT do not include a bypass (although one could be used in the plant loop). The collector requests its design flow rate but it otherwise relies on the larger plant loop for control.

When the PVT thermal collector is controlled to be “on,” in heating mode, and working fluid is flowing, the model calculates the outlet temperature based on the inlet temperature and the collected heat using the following equations.

$$Q_{therm}=A_{surf}\cdot f_{activ}\cdot G_T\cdot {\eta }_{thermal}$$

where:

$Q_{therm}$ is the thermal energy collected (W)

$A_{surf}$ is the net area of the surface (m${}^2$)

$f_{activ}$ is the fraction of surface area with active PV/T collector

${\eta }_{thermal}$ is the thermal conversion efficiency.

$$T_{out}=T_{in}+\frac{Q_{therm}}{\dot{m}{c}_p}$$

where:

$T_{out}$ is the temperature of the working fluid leaving the PV/T

$T_{in}$ is the temperature of the working fluid entering the PV/T

$\dot{m}$ is the entire mass flow of the working fluid through the PV/T

$c_p$ is the specific heat of the working fluid.

For air-based systems, the value of $T_{out}$ is then compared to the temperature setpoint on the outlet node. If $T_{out}$ exceeds the desired outlet temperature, $T_{set,out}$, then a bypass fraction is calculated to model a modulating bypass damper using:

$$f_{bypass}=\frac{\left(T_{set,out}-T_{out}\right)}{\left(T_{in}-T_{out}\right)}$$

When the PVT thermal collector is controlled to be “on,” in cooling mode, and working fluid is flowing, the model calculates the outlet temperature based on the inlet temperature and the heat radiated and convected to the ambient using a heat balance on the outside face of the collector:

$$\dot{m}{c}_p\left(T_{in}-T_{out}\right)={\dot{Q}}_{LWR}+{\dot{Q}}_{conv}$$

where:

${\dot{Q}}_{LWR}$ is the net rate of long wavelength (thermal) radiation exchange with the air, night sky, and ground. See the section “External Longwave Radiation” in the Outside Surface Heat Balance, for full discussion of how this modeled in EnergyPlus using linearized radiation coefficients.

${\dot{Q}}_{conv}$ is the net rate of convective flux exchange with outdoor air. See the section “Exterior/External Convection” in the Outside Surface Heat Balance, for full discussion of how this modeled in EnergyPlus. The surface roughness is assumed to be “very smooth.”

The simple model assumes that the effective collector temperature, $T_{col}$ , is the average of the working fluid inlet and outlet temperatures so that we can make the following substitution:

$$T_{out}=2T_{col}-T_{in}$$

Substituting and solving for $T_{col}$ we obtain the following model for collector temperatures during a (possible) cooling process :

$$T_{col}=\frac{2\dot{m}{c}_p{T}_{in}+A_{surf}{f}_{activ}\left(h_{r,gnd}{T}_{gnd}+h_{r,sky}{T}_{sky}+h_{r,air}{T}_{air}+h_{c,ext}{T}_{air}\right)}{2\dot{m}{c}_p+A_{surf}{f}_{activ}\left(h_{r,gnd}+h_{r,sky}+h_{r,air}+h_{c,ext}\right)}$$

Then the outlet temperature can be calculated and heat losses determined. However, the model allows only sensible cooling of the air stream and limits the outlet temperature to not go below the dewpoint temperature of the inlet.

PVT collectors have a design volume flow rate for the working fluid that is autosizable. For air-based systems used as pre-conditioners, the volume flow rate is sized to meet the maximum outdoor air flow rate. For water-based systems on the supply side of a plant loop, each of the PVT collectors are sized to the overall loop flow rate. For water-based systems on the demand side of a plant loop, the collectors are sized using a rule-of-thumb for typical flow rates per unit of collector area. This rule-of-thumb is based on a constant factor of 1.905x10${}^{-5}$ m${}^3$/s-m${}^2$ that was developed by analyzing SRCC data set for conventional solar collectors (see data set SolarCollectors.idf) and averaging the ratio for all 171 different collectors.

References[LINK]

Charalambous, P.G., Maidment, G.G., Kalagirou, S.A., and Yiakoumetti, K., Photovoltaic thermal (PV/T) collectors: A review. Applied Thermal Engineering 27 (2007) 275-286.

Zondag, H.A. 2008. Flat-plate PV-Thermal collectors and systems: A review. Renewable and Sustainable Energy Reviews 12 (2008) 891-959.

Unglazed Transpired Solar Collectors[LINK]

The input object SolarCollector:UnglazedTranspired provides a model for transpired collectors that are perhaps one of the most efficient ways to collect solar energy with demonstrated instantaneous efficiencies of over 90% and average efficiencies of over 70%. They are used for preheating outdoor air needed for ventilation and processes such as crop drying.

In EnergyPlus, an unglazed transpired solar collector (UTSC) is modeled as a special component attached to the outside face of a heat transfer surface that is also connected to the outdoor air path. A UTSC affects both the thermal envelope and the HVAC air system. From the air system’s point of view, a UTSC is heat exchanger and the modeling needs to determine how much the device raises the temperature of the outdoor air. From the thermal envelope’s point of view, the presence of the collector on the outside of the surface modifies the conditions experienced by the underlying heat transfer surfaces. EnergyPlus models building performance throughout the year and the UTSC will often be “off” in terms of forced airflow, but the collector is still present. When the UTSC is “on” there is suction airflow that is assumed to be uniform across the face. When the UTSC is “off” the collector acts as a radiation and convection baffle situated between the exterior environment and the outside face of the underlying heat transfer surface. We distinguish these two modes of operation as active or passive and model the UTSC component differently depending on which of these modes it is in.

Heat Exchanger Effectiveness[LINK]

The perforated absorber plate is treated as a heat exchanger and modeled using a traditional effectiveness formulation. The heat exchanger effectiveness, ${\epsilon }_{HX}$ , is determined from correlations derived from small-scale experiments. Two correlations available in the literature are implemented in EnergyPlus. The first is based on research by Kutscher at the National Renewable Energy Laboratory. The second is based on the research by Van Decker, Hollands, and Brunger at the University of Waterloo. Because both correlations are considered valid, the choice of which correlation to use is left to the user.

Kutscher Correlation[LINK]

Kutscher’s (1994) correlation encompasses surface convection between the collector and the incoming outdoor air stream that occurs on the front face, in the holes, and along the back face of the collector. The correlation uses a Reynolds number based on the hole diameter as a length scale and the mean velocity of air as it passes through the holes as the velocity scale:

$$R{e}_D=\frac{V_{h}D}{\nu }$$

where:

$V_h$ is the velocity through the holes (m/s)

$D$ is the hole diameter (m)

$\nu$ is the kinematic viscosity of air (m${}^2$/s).

The correlation is a function of Reynolds number, hole geometry, the free stream air velocity, and velocity through the holes:

$$N{u}_D=2.75\left[{\left(\frac{P}{D}\right)}^{-1.2}{\mathrm{Re}}_D^{0.43}+0.011\sigma {\mathrm{Re}}_D{\left(\frac{U_{\infty }}{V_h}\right)}^{0.48}\right]$$

where:

$P$ is the pitch, or distance between holes (m)

$D$ is the diameter of the holes (m)

$\sigma$ is the porosity, or area fraction of the holes

$V_h$ is the mean velocity of air passing through the holes (m/s)

$U_{\infty }$ is the free stream velocity or the local wind speed (m/s).

The Nusselt number is formulated as:

$$N{u}_D=\frac{UD}{k}$$

where:

$U$ is the overall heat transfer coefficient based on log mean temperature difference (W/m${}^2$-K)

$k$ is the thermal conductivity of air (W/m-K).

The heat exchanger effectiveness is:

$${\epsilon }_{HX}=1-e^{\left[-\frac{UA}{\dot{m}{c}_p}\right]}$$

Kutscher’s relation was formulated for triangular hole layout, but based on Van Decker et al. (2001) we allow using the correlation for square hole layout and scale $P$ by a factor of 1.6.

Van Decker, Hollands, and Brunger Correlation[LINK]

Van Decker et. al. extended Kutscher’s measurements to include a wider range of collector parameters including plate thickness, pitch, suction velocities, and square hole patterns. Their model formulation differs from Kutscher’s in that the model was built up from separate effectiveness models for the front, back, and holes of the collector. Their published correlation is:

$$\begin{array}{cc}{\epsilon }_{HX}=& \left[1-\left(1+{\mathrm{Re}}_{s}Max{\left(1.733{\mathrm{Re}}_w^{-\text{1}/12\text{2}},0.02136\right)}^{-1}\right)\right]\times \left[1-{\left(1+0.2273{\mathrm{Re}}_b^{\text{1}/12\text{2}}\right)}^{-1}\right]\\ & \times e^{\left(-0.01895\frac{P}{D}-\frac{20.62}{{\mathrm{Re}}_D}\frac{t}{D}\right)}\end{array}$$

where:

$$R{e}_s=\frac{V_{s}P}{v}$$

$$R{e}_w=\frac{U_{\infty }P}{v}$$

$$R{e}_b=\frac{V_{h}P}{v}$$

$V_s$ is the average suction velocity across the front face of the collector (m/s)

$t$ is the collector plate thickness (m).

Heat Exchanger Leaving Temperature[LINK]

Using either of the correlations above allows determining the heat exchanger effectiveness from known values. By definition the heat exchanger effectiveness is also:

$${\epsilon }_{HX}=\frac{T_{a,HX}-T_{amb}}{T_{s,coll}-T_{amb}}$$

where:

$T_{a,HX}$ is the temperature of the air leaving the collector and entering the plenum (${}^{\circ }$C)

$T_{s,coll}$ is the temperature of the collector’s absorber plate (${}^{\circ }$C)

$T_{amb}$ is the temperature of the ambient outdoor air (${}^{\circ }$C).

By rewriting Equation [eq:SolarCollectorHXEffect716] to solve for $T_{a,HX}$, we see that the temperature of the heated outdoor air entering the plenum can be determined once the collector surface temperature is known:

$$T_{a,HX}={\epsilon }_{HX}{T}_{s,coll}+\left(1-{\epsilon }_{HX}\right)T_{amb}$$

Collector Heat Balance[LINK]

The collector is assumed to be sufficiently thin and high-conductivity so that it can be modeled using a single temperature (for both sides and along its area). This temperature $T_{s,coll}$ is determined by formulating a heat balance on a control volume that just encapsulates the collector surface. The heat balances are formulated separately for active and passive modes and are diagrammed in the following figure.

Observe that for the passive case, we do not use the heat exchanger relations to directly model the interaction of ventilating air with the collector. This is because these relations are considered to not apply when the UTSC is in passive mode. They were developed for uni-directional flow (rather than the balanced-in-and-out flow expected from natural forces) and for specific ranges of suction face velocity. Therefore, this heat transfer mechanism is handled using classical surface convection models (as if the collector was not perforated). (Air exchanges are modeled as ventilation in the plenum air heat balance but do not interact with the hole edges in the collector surface.)

When the UTSC is active, the heat balance on the collector surface control volume is:

$$q_{\alpha sol}^{\prime \prime }+q_{LWR,Env}^{\prime \prime }+q_{conv,wind}^{\prime \prime }-q_{HX}^{\prime \prime }+q_{LWR,plen}^{\prime \prime }+q_{source}^{\prime \prime }=0$$

where:

$q_{\alpha sol}^{\prime \prime }$ is absorbed direct and diffuse solar (short wavelength) radiation heat flux.

$q_{LWR,Env}^{\prime \prime }$ is net long wavelength (thermal) radiation flux exchange with the air and surroundings.

$q_{conv,wind}^{\prime \prime }$ is surface convection flux exchange with outdoor air under high wind and rain conditions. Note that this term is usually assumed to be zero in UTSC model development but we add the term to allow for deteriorated performance of the UTSC under poor conditions.

$q_{HX}^{\prime \prime }$ is heat exchanger flux from collector to incoming outdoor air.

$q_{LWR,plen}^{\prime \prime }$ is net long wavelength (thermal) radiation flux exchange with the outside face of the underlying surface(s).

$q_{source}^{\prime \prime }$ is a source/sink term that accounts for energy exported out of the control volume when the collector’s absorber plate is a hybrid device such as a photovoltaic panel.

While the heat balance on the passive collector surface control volume is:

$$q_{\alpha sol}^{\prime \prime }+q_{LWR,Env}^{\prime \prime }+q_{conv,Env}^{\prime \prime }+q_{LWR,plen}^{\prime \prime }+q_{conv,plen}^{\prime \prime }+q_{source}^{\prime \prime }=0$$

where:

$q_{conv,Env}^{\prime \prime }$ is the surface convection flux exchange with outdoor air

$q_{conv,plen}^{\prime \prime }$ is the surface convection flux exchange with plenum air.

All terms are positive for net flux to the collector except the heat exchanger term, which is taken to be positive in the direction from collector to incoming air stream. Each of these heat balance components is introduced briefly below.

External SW Radiation[LINK]

$q_{\alpha sol}^{\prime \prime }$ is calculated using procedures presented elsewhere in this manual and includes both direct and diffuse incident solar radiation absorbed by the surface face. This is influenced by location, surface facing angle and tilt, shading surfaces, surface face material properties, weather conditions, etc.

External LW Radiation[LINK]

$q_{LWR,Env}^{\prime \prime }$ is a standard radiation exchange formulation between the surface, the sky, the ground, and the atmosphere. The radiation heat flux is calculated from the surface absorptivity, surface temperature, sky, air, and ground temperatures, and sky and ground view factors. Radiation is modeled using linearized coefficients.

External Convection[LINK]

$q_{conv,Env}^{\prime \prime }$ is modeled using the classical formulation:

$$q_{conv}^{\prime \prime }=h_{co}(T_{air}-T_o)$$

where h${}_{co}$ is the convection coefficient. This coefficient will differ depending on whether or not the UTSC is active or passive. When the UTSC is passive, h${}_{co}$ is treated in the same way as an outside face with ExteriorEnvironment conditions. When the UTSC is active, the special suction airflow situation of a transpired collector during operation means that h${}_{co}$ is often zero because the suction situation can eliminate mass transport away from the collector. However when the winds are high, the strong turbulence and highly varying pressures can cause the suction flow situation to breakdown. Therefore, we include the $q_{conv,wind}^{\prime \prime }$ term in the heat balance and use a special coefficient $h_{c,wind}$ to model this lost heat transfer. In addition, when it is raining outside, we assume the collector gets wet and model the enhanced surface heat transfer using a large value for $h_{c,wind}$.

Heat Exchanger[LINK]

$q_{HX}^{\prime \prime }$ is modeled using the classical formulation:

$$q_{HX}^{\prime \prime }=\frac{\dot{m}{c}_p(T_{a,HX}-T_{amb})}{A}$$

where $T_{a,HX}$ is determined using correlations described above. When the UTSC is active, the air mass flow is determined from the operation of the outdoor air mixer component. When the UTSC is off, this term is zero.

Plenum LW Radation[LINK]

$q_{LWR,plen}^{\prime \prime }$ is a standard radiation exchange formulation between the collector surface and the underlying heat transfer surface located across the plenum. Radiation is modeled using linearized coefficients.

Plenum Convection[LINK]

$q_{conv,plen}^{\prime \prime }$ is modeled using the classical formulation:

$$q_{conv}^{\prime \prime }=h_{cp}(T_{air}-T_o)$$

where h${}_{cp}$ is the convection coefficient. This coefficient is taken as zero when the UTSC is operating because of the suction airflow situation. When the UTSC is off, the value for h${}_{cp}$ is obtained from correlations used for window gaps from ISO (2003) standard 15099.

Substituting these models into Equation [eq:SolarCollectorHeatBalance718] and solving for $T_{s,coll}$ yields the following equation when the UTSC is active (“on”):

$$T_{s,coll}=\frac{\left(I_s\alpha +h_{r,atm}{T}_{amb}+h_{r,sky}{T}_{sky}+h_{r,gnd}{T}_{amb}+h_{r,plen}{T}_{so}+h_{c,wind}{T}_{amb}+\frac{\dot{m}{c}_p}{A}{T}_{amb}-\frac{\dot{m}{c}_p}{A}\left(1-{\epsilon }_{HX}\right)T_{amb}+{q^{\prime \prime }}_{source}\right)}{\left(h_{r,atm}+h_{r,sky}+h_{r,gnd}+h_{r,plen}+h_{c,wind}+\frac{\dot{m}{c}_p}{A}{\epsilon }_{HX}\right)}$$

and substituting into Equation [eq:SolarCollectorHeatBalance718] yields the following equation when the UTSC is passive (“off”):

$$T_{s,coll}=\frac{\left(I_s\alpha +h_{co}{T}_{amb}+h_{r,atm}{T}_{amb}+h_{r,sky}{T}_{sky}+h_{r,gnd}{T}_{amb}+h_{r,plen}{T}_{so}+h_{c,plen}{T}_{a,plen}+{q^{\prime \prime }}_{source}\right)}{\left(h_{co}+h_{r,air}+h_{r,sky}+h_{r,gnd}+h_{r,plen}+h_{c,plen}\right)}$$

where:

$I_s$ is the incident solar radiation of all types (W/m${}^2$)

$\alpha$ is the solar absorptivity of the collector

$h_{r,atm}$ is the linearized radiation coefficient for the surrounding atmosphere (W/m${}^2$-K)

$T_{amb}$ is the outdoor drybulb from the weather data, also assumed for ground surface (${}^{\circ }$C)

$h_{r,sky}$ is the linearized radiation coefficient for the sky (W/m${}^2$-K)

$T_{sky}$ is the effective sky temperature (${}^{\circ }$C)

$h_{r,gnd}$ is the linearized radiation coefficient for the ground (W/m${}^2$-K)

$h_{r,plen}$ is the linearized radiation coefficient for the underlying surface (W/m${}^2$-K)

$h_{c,wind}$ is the convection coefficient for the outdoor environment when the UTSC is active and winds are high or it is raining (W/m${}^2$-K)

$T_{so}$ is the temperature of the outside face of the underlying heat transfer surface (${}^{\circ }$C)

$\dot{m}$ is the air mass flow rate when in active mode (kg/s)

$c_p$ is the specific heat of air at constant pressure (J/kg-K)

$A$ is the area of the collector (m${}^2$)

$h_{co}$ is the convection coefficient for the outdoor environment (W/m${}^2$-K)

$h_{c,plen}$ is the convection coefficient for the surfaces facing the plenum (W/m${}^2$-K)

$T_{a,plen}$ is the air drybulb for air in the plenum and entering the outdoor air system (${}^{\circ }$C).

Plenum Heat Balance[LINK]

The plenum is the volume of air located between the collector and the underlying heat transfer surface. The plenum air is modeled as well-mixed. The uniform temperature of the plenum air, $T_{a,plen}$, is determined by formulating a heat balance on a control volume of air as diagrammed below.

Note that we have formulated the control volumes with slight differences for the active and passive cases. For the active case, the suction air situation and heat exchanger effectiveness formulations dictate that the collector surface control volume encompass part of the air adjacent to both the front and back surfaces of the collector. However for the passive case, the collector surface control volume has no air in it and the plenum air control volume extends all the way to the surface of the collector.

When the UTSC is active, the heat balance on the plenum air control volume is:

$${\dot{Q}}_{air}+{\dot{Q}}_{co}=0$$

where:

${\dot{Q}}_{air}$ is the net rate of energy added by suction air convecting through the control volume

${\dot{Q}}_{co}$ is the net rate of energy added by surface convection heat transfer with the underlying surface.

When the UTSC is passive, the heat balance on the plenum air control volume is:

$${\dot{Q}}_{vent}+{\dot{Q}}_{co}+{\dot{Q}}_{c,coll}=0$$

where:

${\dot{Q}}_{vent}$ is the net rate of energy added from infiltration – where outdoor ambient air exchanges with the plenum air

${\dot{Q}}_{c,coll}$ is the net rate of energy added by surface convection heat transfer with the collector.

Substituting into Equation [eq:SolarCollectorHeatBalancePlenumAir722] and solving for $T_{a,plen}$ yields the following equation for when the UTSC is active:

$$T_{a,plen}=\frac{\left(\dot{m}{c}_p{T}_{a,HX}+h_{c,plen}A{T}_{so}\right)}{\left(\dot{m}{c}_p+h_{c,plen}A\right)}$$

And substituting into Equation [eq:SolarCollectorHeatBalancePlenumAir722] yields the following equation when the UTSC is passive:

$$T_{a,plen}=\frac{\left(h_{c,plen}A{T}_{so}+{\dot{m}}_{vent}{c}_p{T}_{amb}+h_{c,plen}A{T}_{s,coll}\right)}{\left(h_{c,plen}A+{\dot{m}}_{vent}{c}_p+h_{c,plen}A\right)}$$

where ${\dot{m}}_{vent}$ is the air mass flow from natural forces (kg/s).

The literature on UTSC does not appear to address the passive mode of operation and no models for ${\dot{m}}_{vent}$ have been identified. Nevertheless, natural buoyancy and wind forces are expected to drive air exchange between the plenum and ambient and some method of modeling ${\dot{m}}_{vent}$ is needed. Reasoning that the configuration is similar to single-side natural ventilation, we elect to use correlations for natural ventilation presented in Chapter 26 of the ASHRAE Handbook of Fundamentals (2001).

$${\dot{m}}_{vent}=\rho {--\dot{V}}_{tot}$$

where:

$\rho$ is the density of air (kg/m${}^3$)

${--\dot{V}}_{tot}={--\dot{V}}_{wind}+{--\dot{V}}_{thermal}$ is the total volumetric flow rate of air ventilating in and out of the plenum (m${}^3$/s)

$${--\dot{V}}_{wind}=C_v{A}_{in}{U}_{\infty }$$

$${--\dot{V}}_{thermal}=C_D{A}_{in}\sqrt{2g\Delta H_{NPL}\left(T_{a,plen}-T_{amb}\right)/T_{a,plen}}(if{T}_{a,plen}>T_{amb})$$

$${--\dot{V}}_{thermal}=C_D{A}_{in}\sqrt{2g\Delta H_{NPL}\left(T_{amb}-T_{a,plen}\right)/T_{amb}}(if{T}_{amb}>T_{a,plen}andUTSCvertical)$$

where:

$C_v$ is the effectiveness of the openings that depends on opening geometry and the orientation with respect to the wind. ASHRAE HoF (2001) indicates values ranging from 0.25 to 0.6. In the UTSC model, this value is available for user input and defaulted to 0.25.

$C_D$ is the discharge coefficient for the opening and depends on opening geometry. In the UTSC model, this value is available for user input and defaulted to 0.65.

Mass continuity arguments lead to modeling the area of the openings as one half of the total area of the holes, so we have:

$$A_{in}=\frac{A\sigma }{2}$$

$g$ is the gravitational constant taken as 9.81 (m/s${}^2$).

$\Delta H_{NPL}$ is the height from midpoint of lower opening to the Neutral Pressure Level. This is taken as one-fourth the overall height of the UTSC if it is mounted vertically. For tilted collectors, the nominal height is modified by the sine of the tilt. If the UTSC is mounted horizontally (e.g. on the roof) then the $\Delta H_{NPL}$ is taken as the gap thickness of the plenum.

If the UTSC is horizontal and $T_{amb}>T_{a,plen}$, then ${--\dot{V}}_{thermal}=0$ because this is a stable situation.

Underlying Heat Transfer Surface[LINK]

The UTSC is applied to the outside of a heat transfer surface. This surface is modeled using the usual EnergyPlus methods for handling heat capacity and transients – typically the CTF method. These native EnergyPlus Heat Balance routines are used to calculate $T_{so}$ . The UTSC model is coupled to the underlying surface using the OtherSideConditionsModel mechanism. The UTSC model provides values for $h_{r,plen}$, $T_{s,coll}$, $h_{c,plen}$, and $T_{a,plen}$ for use with the Heat Balance Model calculations for the outside face of the underlying surface (described elsewhere in this manual).

Solar and Shading Calculations[LINK]

The transpired collector object uses a standard EnergyPlus surface in order to take advantage of the detailed solar and shading calculations. Solar radiation incident on the surface includes beam and diffuse radiation, as well as radiation reflected from the ground and adjacent surfaces. Shading of the collector by other surfaces, such as nearby buildings or trees, is also taken into account.

Local Wind Speed Calculations[LINK]

The outdoor wind speed affects terms used in modeling UTSC components. The wind speed in the weather file is assumed to be measured at a meteorological station located in an open field at a height of 10 m. To adjust for different terrain at the building site and differences in the height of building surfaces, the local wind speed is calculated for each surface.

The wind speed is modified from the measured meteorological wind speed by the equation (ASHRAE 2001):

$$U_{\infty }=V_{met}{\left(\frac{{\delta }_{met}}{z_{met}}\right)}^{a_{met}}{\left(\frac{z}{\delta }\right)}^a$$

where:

$z$ is the height of the centroid of the UTSC

$z_{met}$ is the height of the standard meteorological wind speed measurement

$a$ and $T_{va}$ are terrain-dependent coefficients. $T_{va}$ is the boundary layer thickness for the given terrain type. The values of a and $T_{va}$ are shown in the Table 1.

The UTSC can be defined such that it has multiple underlying heat transfer surfaces. The centroid heights for each surface are area-weighted to determine the average height for use in the local wind calculation.

Convection Coefficients[LINK]

UTSC modeling requires calculating up to three different coefficients for surface convection heat transfer. These coefficients are defined in the classic way by:

$$h_c=\frac{T_{air}-T_{surf}}{{q^{\prime \prime }}_{conv}}$$

First, $h_{co}$ is the convection coefficient for the collector surface facing the outdoors when the UTSC is passive. It is modeled in exactly the same way as elsewhere in EnergyPlus and will depend on the user setting for Outside Convection Algorithm – Outside Surface Heat Balance entry elsewhere in this document.

Second, $h_{c,plen}$ is the convection coefficient for surfaces facing the plenum. This coefficient is applied to just the underlying surface’s convection when the UTSC is active and to both the collector and the underlying surface when the UTSC is passive. When the UTSC is active, we use the convection correlation for forced air developed by McAdams (1954) as published by ASHRAE HoF (2001):

$$h_{c,plen}=5.62+3.9V_p$$

where $V_p$ is the mean velocity in the plenum determined from:

$$V_p=\frac{\dot{m}}{2\rho A_p}$$

and where $A_p$ is the effective cross section area of the plenum perpendicular to the primary flow direction. When the UTSC is passive, we model the convection in the same way used in EnergyPlus to model air gaps in windows. These correlations vary by Rayleigh number and surface tilt and are based on the work of various research including Hollands et. al., Elsherbiny et. al., Wright, and Arnold. The formulations are documented in ISO (2003) standard 15099. For the UTSC implementation, the routines were adapted from Subroutine NusseltNumber in WindowManager.f90 (by F. Winkelmann), which itself was derived from Window5 subroutine “nusselt”.

Third, $h_{c,wind}$ is the convection coefficient used to degrade the UTSC performance under environmental conditions with high wind or rain. If the weather file indicates it is raining, then we set $h_{c,wind}$ = 1000.0 which has the effect of making the collector the same temperature as ambient air. The heat exchanger effectiveness correlations described above account for a moderate amount of wind, but the correlations appear confined to the range 0 to 5.0 m/s. Therefore we set $h_{c,wind}$ equal to zero if $U_{\infty }$ is < = 5.0 m/s. If $U_{\infty }$ is > 5.0 m/s, then we use the McAdams correlation but with a reduced velocity magnitude:

$$h_{c,wind}=5.62+3.9(U_{\infty }-5.0)$$

Radiation Coefficients[LINK]

UTSC modeling requires calculating up to four different linearized coefficients for radiation heat transfer. Whereas radiation calculations usually use temperature raised to the fourth power, this greatly complicates solving heat balance equations for a single temperature. Linearized radiation coefficients have the same units and are used in the same manner as surface convection coefficients and introduce very little error for the temperature levels involved.

The radiation coefficient, $h_{r,plen}$, is used to model thermal radiation between the collector surface and the outside face of the underlying heat transfer surface. We assume a view factor of unity. It is calculated using:

$$h_{r,plen}={\sigma }_{SB}{e}_{coll}{e}_{so}\frac{\left(T_{s,coll}^4-T_{so}^4\right)}{\left(T_{s,coll}-T_{so}\right)}$$

where all temperatures are converted to Kelvin and:

${\sigma }_{SB}$ is the Stefan-Boltzmann constant

$e_{coll}$ is the longwave thermal emittance of the collector

$e_{so}$ is the longwave thermal emittance of the underlying heat transfer surface.

The three other coefficients, $h_{r,atm}$, $h_{r,sky}$, and $h_{r,gnd}$ are used elsewhere in EnergyPlus for the outside face surface heat balance and are calculated in the same manner as equation for UTSC collectors. [This is accomplished by calling subroutine InitExteriorConvectionCoeffs in the file HeatBalanceConvectionCoeffs.f90.]

Bypass Control[LINK]

The UTSC is assumed to be arranged so that a bypass damper controls whether or not air is drawn directly from the outdoors or through the UTSC. The control decision is based on whether or not it will be beneficial to heat the outdoor air. There are multiple levels of control including an availability schedule, whether or not the outdoor air is cooler than the mixed air setpoint, or whether or not the zone air temperature is lower than a so-called free heating setpoint.

Sizing Warnings[LINK]

Although the design of the transpired collector is left to the user, the program issues warnings when the suction airflow velocity falls outside the range 0.003 to 0.08 m/s.

Overall Efficiency[LINK]

The overall thermal efficiency of the UTSC is a useful output report and is defined as the ratio of the useful heat gain of the entire system versus the total incident solar radiation on the gross surface area of the collector.

$$\eta =\frac{\left(\dot{Q}/A\right)}{I_{sc}}=\frac{\dot{m}{c}_p\left(T_{a,plen}-T_{amb}\right)}{I_{sc}A}$$

where:

$\dot{Q}$ is useful heat gain

$I_{sc}$ is total incident solar radiation.

Note that the efficiency $\eta$ is only defined for $I_{solar}>0$. This efficiency includes heat recovered from the underlying wall and can exceed 1.0.

Collector Efficiency[LINK]

The thermal efficiency of the collector is a useful output report and is defined as the ratio of the useful heat gain of the collector fluid versus the total incident solar radiation on the gross surface area of the collector.

$$\eta =\frac{\dot{m}{c}_p\left(T_{a,HX}-T_{amb}\right)}{I_{sc}A}$$

Note that the efficiency $\eta$ is only defined for $I_{solar}>0$.

References[LINK]

Kutscher, C.F. 1994. Heat exchange effectiveness and pressure drop for air flow through perforated plates with and without crosswind. Journal of Heat Transfer. May 1994, Vol. 116, p. 391. American Society of Mechanical Engineers.

Van Decker, G.W.E., K.G.T. Hollands, and A.P. Brunger. 2001. Heat-exchange relations for unglazed transpired solar collectors with circular holes on a square of triangular pitch. Solar Energy. Vol. 71, No. 1. pp 33-45, 2001.

ISO. 2003. ISO 15099:2003. Thermal performance of windows, doors, and shading devices – Detailed calculations. International Organization for Standardization.

Air-Based Building Integrated Parallel Flow PV/Thermal Collector[LINK]

This section discusses the implementation of an air-based parallel flow building-integrated PV/Thermal collector model. These collectors are directly integrated in the building façade or roof. The photovoltaic panels become the cladding material. Air is blown in the cavity behind the photovoltaic panels to pick up heat from the photovoltaic panels. This heat usually helps to reduce ventilation air heating requirements during the heating season. It also helps to marginally improve the PV efficiency by keeping the cells cool.

Energy Balance[LINK]

The input object SolarCollector:FlatPlate:PhotovoltaicThermal provides models for flat plate photovoltaic/thermal collectors. To use the BIPVT collector model, use the SolarCollectorPerformance:PhotovoltaicThermal:BIPVT object and specify its name under the parameter ‘Photovoltaic-Thermal Model Performance Name’ of the SolarCollector:FlatPlate:PhotovoltaicThermal object. The heat balances for the various layers of the BIPVT collector discussed here is taken from work by Delisle and Kummert (2014), but the solution methodology is different (discussed in a later section). The thermal resistance circuit representation is shown in Figure 7. The equations [eq:BIPVT-pv-outer-surf] through [eq:BIPVT-cavity-air] describe energy balances on the various layers of interest. These are the PV glazing outer layer, the PV cells, the PV backing outer surface, and the air cavity. The energy balances are for a length (along flow direction) of collector $dx$.

PV glazing outer surface: $$0=h_{conv,t}(T_{PVg}-T_a)+h_{rad,t}(T_{PVg}-T_{sur})+\left(\frac{T_{PVg}-T_{PV}}{R_{PVg-PV}}\right)$$

PV cells: $$S+\left(\frac{T_{PVg}-T_{PV}}{R_{PVg-PV}}\right)=\left(\frac{T_{PV}-T_1}{R_{PV-1}}\right)$$

PV backing outer surface: $$\left(\frac{T_{PV}-T_1}{R_{PV-1}}\right)+S_1=h_{conv,f1}(T_1-T_f)+h_{rad,1-2}(T_1-T_2)$$

Fluid (air) in cavity: $$\dot{m}{C}_p\frac{d{T}_f}{dx}=[h_{conv,f1}(T_1-T_f)+h_{conv,f2}(T_2-T_f)]\cdot W$$

where:

$$S=G_{PV}\cdot \tau \alpha (0{)}_{PV}\cdot F_{cell}\cdot F_{module}-{\eta }_{PV}\cdot G_{PV}\cdot F_{module}$$

$$S_1=G_{bs}\cdot \tau \alpha (0{)}_{bs}\cdot (1-F_{cell})\cdot F_{module}+\tau \alpha (0{)}_{clad}\cdot G_{clad}\cdot (1-F_{module})$$

$$G_{PV}=G_{beam}\cdot IA{M}_{beam,PV}+G_{sky}\cdot IA{M}_{sky,PV}+G_{ground}\cdot IA{M}_{ground,PV}$$

$$G_{bs}=G_{beam}\cdot IA{M}_{beam,bs}+G_{sky}\cdot IA{M}_{sky,bs}+G_{ground}\cdot IA{M}_{ground,bs}$$

$$G_{clad}=G_{beam}\cdot IA{M}_{beam,clad}+G_{sky}\cdot IA{M}_{sky,clad}+G_{ground}\cdot IA{M}_{ground,clad}$$

Note: For simplicity, incident angle modifier of the PV backing material, cladding, and of the PV cell are assumed to be the same here, so $G_{bs}=G_{clad}=G_{PV}$. This is a generally good assumption, as can be calculated using work by Beckman et al. (1977) for a flat black surface, presented in Duffie and Beckman (2006), and comparing to the PV/glass interface model using Snell’s, Fresnel’s and Bouguer’s law (described in the next section). The results for both are very similar.

$$\begin{array}{cc}{T}_{sur}=& \\ & {[{(T_a+273.15)}^4\cdot 0.5\cdot (1-\cos \left(\beta \right))+{(T_{sky}+273.15)}^4\cdot 0.5\cdot (1+\cos \left(\beta \right))]}^{1/4}-273.15\end{array}$$

$T_{f,i}$ is specific heat of the fluid (air) ($J\cdot k{g}^{-1}{K}^{-1}$)

$F_{cell}$ is the area fraction of cells on PV modules (i.e. percentage of PV module area that is covered with PV cells) (-)

$F_{module}$ is the area fraction of PV modules on wall/roof surface (i.e. percentage of roof or wall area that is covered by PV modules) (-)

$G_{beam}$ is the incident beam radiation on the collector ($W\cdot m^{-2}$)

$G_{sky}$ is the incident sky diffuse radiation on the collector ($W\cdot m^{-2}$)

$G_{ground}$ is the incident ground reflected diffuse radiation on the collector ($W\cdot m^{-2}$)

$h_{conv,t}$ is the convective heat transfer coefficient between the PV glazing and the ambient air ($W\cdot m^{-2}{K}^{-1}$)

$h_{rad,t}$ is the radiative heat transfer coefficient between the PV glazing and the surrounding surfaces ($W\cdot m^{-2}{K}^{-1}$)

$h_{conv,f1}$ is the convective heat transfer coefficient between the air in the cavity and the back of the PV modules ($W\cdot m^{-2}{K}^{-1}$)

$h_{conv,f2}$ is the convective heat transfer coefficient between the air in the cavity and the building wall/roof ($W\cdot m^{-2}{K}^{-1}$)

$h_{rad,1-2}$ is the radiative heat transfer coefficient between the back of the PV modules and the building wall/roof ($W\cdot m^{-2}{K}^{-1}$)

$IA{M}_{a,b}$ is the incident angle modifier for incident radiation a, on surface b (-)

$\dot{m}$ is the flow rate of the fluid (air) ($kg\cdot s^{-1}$)

$R_{PVg-PV}$ is the thermal resistance between the surface PV glazing and the PV cells ($m^{2}K\cdot W^{-1}$)

$R_{PV-1}$ is the thermal resistance between the PV cells and the back of the PV modules ($m^{2}K\cdot W^{-1}$)

$T_{PVg}$ is the temperature of the surface of the PV glazing (${}^{\circ }C$)

$T_{PV}$ is the temperature of the PV cells (${}^{\circ }C$)

$T_a$ is the ambient temperature (${}^{\circ }C$)

$T_{sky}$ is the effective sky temperature (${}^{\circ }C$)

$T_1$ is the temperature of the back of the PV modules (${}^{\circ }C$)

$T_2$ is the temperature of the surface of the building wall/roof (${}^{\circ }C$)

$T_f$ is the temperature of the fluid (air) between the PV module and the building surface (${}^{\circ }C$)

$W$ is the width of the collector ($m$)

$\beta$ is the slope of the collector surface (0${}^{\circ }$ is horizontal, 90${}^{\circ }$ is vertical) (degrees)

${\eta }_{PV}$ is the PV efficiency at current (operating) condition, calculated by Generator:Photovoltaic (-)

$\tau \alpha (0{)}_a$ is the nominal or measured transmittance-absorptance product of surface $a$ at normal incidence angle (-)

Incident angle modifier[LINK]

The incident angle modifiers $IA{M}_{a,b}$ were calculated using Snell’s Fresnel’s and Bouguer’s law, as described in Duffie and Beckman (2006). As mentioned previously, the $IAM$ for the PV cells is assumed to be the very similar to the IAM for other surfaces, so the PV cell IAM is used for all surfaces. This section shows the calculation for the IAM of the PV cells. Equation [eq:BIPVT-IAM] shows the $IAM$ for an incident angle $\theta$. This can be applied to any source of incident solar radiation (beam, sky, ground). Equation [eq:BIPVT-IAM-2] assumed no reflection from the PV cell (hence the subscript $\alpha =1$), which is not true in reality, but simplifies the calculations. The absorption of the cell (or other material) is accounted for in Equations [eq:BIPVT-S] and [eq:BIPVT-S] with the use of the ‘nominal’ or ‘measured’ $\tau \alpha$.

$$IA{M}_{\theta }=\frac{\tau \alpha (\theta {)}_{\alpha =1}}{\tau \alpha (0{)}_{\alpha =1}}$$

where: $$\tau \alpha (\theta {)}_{\alpha =1}=e^{-(K\cdot d/\cos \left({\theta }_r\right))}[1-\frac{1}{2}(\frac{{\sin }^2({\theta }_r-\theta )}{{\sin }^2({\theta }_r+\theta )}+\frac{{\tan }^2({\theta }_r-\theta )}{{\tan }^2({\theta }_r+\theta )})]$$

$${\theta }_r={\sin }^{-1}\left(\sin \right(\theta \left)/n_{glass}\right)$$

$n_{glass}$ is the refractive index of glass (around 1.526) (-)

${\theta }_r$ is the angle of refraction in the glazing (degrees)

$\theta$ is the angle of incidence (degrees)

$K$ is the glazing extinction coefficient ($m^{-1}$)

$d$ is the thickness of the glazing ($m$)

$\tau \alpha (\theta {)}_{\alpha =1}$ and $\tau \alpha (0{)}_{\alpha =1}$ are the calculated transmittance-absorptance products at incidence angle $\theta$ and 0 respectively. For calculating the incident angle modifier, these are calculated assuming that absorptivity of the PV cell is equal to 1.

Conduction heat transfer coefficients[LINK]

Conduction heat transfer coefficients can be calculated as the inverse of their associated resistances as shown below.

$$h_{PVg-PV}=\frac{1}{R_{PVg-PV}}$$

$$h_{PV-1}=\frac{1}{R_{PV-1}}$$

Radiative heat transfer coefficients[LINK]

The radiative heat transfer coefficients are given by:

$$h_{rad,t}={ϵ}_{PVg}\sigma (T_{PVg}^2+T_{sur}^2)(T_{PVg}+T_{sur})$$

$$h_{1,2}=\sigma \frac{(T_1^2+T_2^2)(T_1+T_2)}{\frac{1}{{ϵ}_1}+\frac{1}{{ϵ}_2}-1}$$

where:

$\sigma$ is the Stefan-Boltzmann constant (approx. 5.67E-8 $W\cdot m^{-2}{K}^{-4}$)

${ϵ}_a$ is the long wave emissivity of surface $a$ (-)

$h_{1,2}$ is the radiant heat transfer coefficient between surfaces 1 and 2 ($W\cdot m^{-2}{K}^{-1}$)

$h_{rad,t}$ is the radiant heat transfer coefficient between the top surface of the collector and the surroundings ($W\cdot m^{-2}{K}^{-1}$)

Convective heat transfer coefficients[LINK]

Internal and external heat transfer coefficients have a great impact on the results of unglazed collectors. Therefore, heat transfer coefficients must be carefully chosen, while acknowledging that there will always be a certain level of uncertainty resulting from models like this.

When the flow is turbulent inside the air channel, the heat transfer coefficients $h_{conv,f1}$ and $h_{conv,f2}$ are calculated using the correlations by Candanedo et al. (2011) (Equations [eq:BIPVT-NuTop] and [eq:BIPVT-NuBott]) for the Nusselt Number of the top an bottom surfaces respectively.

$$N{u}_{D_{h,c}}=0.052R{e}^{0.78}P{r}^{0.4}$$

$$N{u}_{D_{h,c}}=1.017R{e}^{0.471}P{r}^{0.4}$$

where:

$Re$ is the Reynold number (-)

$Pr$ is the Prandlt number (-)

$D_{h,c}$ is the hydraulic diameter of the collector channel (-)

The external forced convective heat transfer coefficient is derived from results by Gorman et al (2019). Data from Figure 8 shows graphical results for external heat transfer coefficients. A correlation of the form shown in equation [eq:BIPVT-ext-htc] was used to fit to the data from Gorman et al.

$$h_{conv,ext,forced}=\frac{Y\cdot V^m}{L_{ch}^{1-m}}$$

where:

$L_{ch}$ is the characteristic length (for roof – length along flow direction, for leeward and windward vertical surfaces – hydraulic diameter of surface, vertical sides – length of surfaces along flow direction) ($m$)

$V$ is wind velocity ($m\cdot s^{-1}$)

$Y$ and $m$ are calculated coefficient and exponent

The exponents $m$ derived for the various surfaces are consistent with theory and literature. Laminar flow has a relationship of $Nu\propto {Re}^{0.5}$, fully developed turbulent flow $Nu\propto {Re}^{0.8}$ (Incropera et al, 2007), and fully separated flow $Nu\propto {Re}^{0.67}$(Sparrow et al, 1982). As shown in Table 2, the windward roof has an exponent nearing 0.8 (fully developed, turbulent flow), and the other surfaces have exponents between 0.6 and 0.7 (this could mean that either the flow is separating, or that the flow is not fully developed).

External natural convection is calculated from Incropera et al. (2007) (shown in equation [eq:BIPVT-htc-nat]).

$$h_{conv,ext,nat}=0.15\cdot {Ra}^{1/3}\cdot \frac{k}{D_{h,s}}$$

where:

$Ra$ is Raleigh number (-)

$k$ is thermal conductivity of air ($W\cdot m^{-1}{K}^{-1}$)

$D_{h,s}$ is hydraulic diameter of the heat transfer surface ($m$)

The combined external heat transfer coefficient is calculated with the method described in Incropera et al. (2007) This is shown in equation [eq:BIPVT-htc-conv-ext], with an exponent of 3 which, as explained by Incropera et al. (2007), is a good value for most cases.

$$h_{conv,ext}={({h_{conv,ext,forced}}^3+{h_{conv,ext,nat}}^3)}^{1/3}$$

Solving energy balance equations for positive flow condition[LINK]

To solve the energy balance equations when $\dot{m}>0$ (air is flowing through collector) shown in equations [eq:BIPVT-pv-outer-surf] through [eq:BIPVT-cavity-air], first let: $$A=-\frac{W}{\dot{m}{C}_p}(h_{conv,f1}+h_{conv,f2})$$

$$B=\frac{W}{\dot{m}{C}_p}(h_{conv,f1}\cdot T_1+h_{conv,f2}\cdot T_2)$$

re-arranging equation [eq:BIPVT-cavity-air] yields: $$\frac{d{T}_f}{dx}=A\cdot T_f+B$$

Re-arrange equation [eq:BIPVT-dtdx] and integrate both sides for the entire length of the collector (length=$L$).

$${\int }_{x=0}^{x=L}\frac{d{T}_f}{{A\cdot T}_f+B}={\int }_{x=0}^{x=L}dx$$

For this solution, we assume that $T_1$ and $T_2$ are not a function of distance $x$ (i.e. $T_1$ and $T_2$ are average temperatures over collector length, $\overline{T_1}$ and $\overline{T_2}$). Solving equation [eq:BIPVT-integral] gives equation [eq:BIPVT-Tfl].

$$T_f\left(L\right)=\left(T_f\right(0)+\frac{B}{A})e^{A\cdot L}-\frac{B}{A}$$

If the fluid temperature along another point along the collector is required, replace $L$ with required distance $x$.

The average fluid temperature over the length of the collector is calculated in the following manner:

$$\overline{T_f}=\frac{1}{L}{\int }_{x=0}^{x=L}{T}_f(x)dx=\frac{1}{A\cdot L}(T_f\left(0\right)+\frac{B}{A})(e^{A\cdot L}-1)-\frac{B}{A}$$

We have previously assumed that all other layer’s temperatures are not functions of $x$ (i.e. they are average temperatures over the length of the collector). If we also assume that $T_f\left(0\right)$, $\overline{T_2}$, $S$, and $S_1$ are known, we can solve equations [eq:BIPVT-pv-outer-surf],[eq:BIPVT-pv-cells], and [eq:BIPVT-PV-back]. Equations [eq:BIPVT-pv-outer-surf] through [eq:BIPVT-PV-back] can be re-written in the following form:

$$\begin{array}{c}\left[\begin{array}{ccc}{h}_{conv,t}+h_{rad,t}+h_{PVg-PV}& -h_{PVg-PV}& 0\\ h_{PVg-PV}& -h_{PVg-PV}-h_{PV-1}& h_{PV-1}\\ 0& h_{PV-1}& -h_{PV-1}-h_{conv,f1}-h_{rad,1-2}\end{array}\right]\left[\begin{array}{c}{\overline{T}}_{PVg}\\ {\overline{T}}_{PV}\\ {\overline{T}}_1\end{array}\right]\\ =\left[\begin{array}{c}{h}_{conv,t}{T}_a+h_{rad,t}{T}_{sur}\S {\S }_1-h_{rad,1-2}{\overline{T}}_2-h_{conv,f1}\overline{T_f}\end{array}\right]\end{array}$$

Equation [eq:BIPVT-matrix] can then be solved with matrix inversion.

The PV temperature will be passed directly to the proper Photovoltaics subroutine to calculate the PV efficiency.

The following steps are used to solve the above equations at each time step:

1. Assume (guess) values ${\overline{T}}_f$ and ${\overline{T}}_{PV}$. If not the first time step, use previous time step values.

2. Get current values of all boundary conditions

3. Get PV efficiency ${\eta }_{PV}$ from PV model using ${\overline{T}}_{PV}$ from step 1

4. Update all the coefficients and constants in Equation [eq:BIPVT-matrix]

5. Solve Equation [eq:BIPVT-matrix] (Gives values for ${\overline{T}}_{PVg}$, ${\overline{T}}_{PV}$, ${\overline{T}}_1$)

6. Solve Equation [eq:BIPVT-tfl-avg] (Gives ${\overline{T}}_f$)

7. Iterate steps 2-6 until convergence

8. Calculate $T_f\left(L\right)$ using Equation [eq:BIPVT-Tfl]

9. Calculate $\dot{Q}=\dot{m}{C}_p\left(T_f\right(L)-T_f(0\left)\right)$

Solving energy balance equations for zero flow condition[LINK]

A similar approach is taken when the collector is stagnating. Calculations of the temperatures in the collector is required to accurately assess the PV efficiency even though no thermal output is generated during stagnation.

In stagnation conditions, Equations [eq:BIPVT-pv-outer-surf]-[eq:BIPVT-cavity-air] with $\dot{m}=0$ are summed to give:

$$T_{f,avg}=-\frac{D}{E}$$

where:

$$C=S+S_1+h_{conv,t}\left({T_a-T}_{PVg}\right)+h_{rad,t}(T_{sur}-T_{PVg})+h_{rad,1-2}(T_2-T_1)$$

$$D=C+h_{conv,f2}\cdot T_2$$

$$E=-h_{conv,f2}$$

We use the correlation by Hollands et al (1976) for the entire cavity heat transfer coefficient (i.e. heat transfer coefficient between the plates) as shown in equation [eq:BIPVT-h-cavity]. Since the correlation by Hollands et al was for the entire cavity, and if we assume both that the internal heat transfer coefficients $h_{conv,f2}$ and $h_{conv,f1}$ are equal, we can calculate $h_{conv,f2}$ using equation [eq:BIPVT-h-cavity-2].

$$h_{cavity}=\frac{k_{air}}{d}(1+1.44\cdot [1-\frac{1708\cdot {\left(\sin 1.8\beta \right)}^{1.6})}{R{a}_d\cos \beta }]{[1-\frac{1708}{R{a}_d\cos \beta }]}^{\ast }+{[{\left(\frac{R{a}_d\cos \beta }{5830}\right)}^{1/3}-1]}^{\ast })$$

where:

$$Ra=\frac{g\left(\frac{1}{T_f+273.15}\right)|T_2-T_1|d^3}{\nu \alpha }$$

$[{]}^{\ast }$ denotes that if the value in brackets is must be positive, else be set to zero.

$d$ is the thickness of the channel ($m$)

$\upsilon$ is the thermal conductivity of air in the cavity ($W\cdot m^{-1}{K}^{-1}$)

$g$ is gravity ($m\cdot s^2$)

$\nu$ is the kinematic viscosity of air ($m^2{s}^{-1}$)

$\alpha$ is the thermal diffusivity of air ($m^2{s}^{-1}$)

$$\frac{1}{h_{cavity}}=\frac{1}{h_{conv,f1}}+\frac{1}{h_{conv,f2}}=\frac{2}{h_{conv,f2}}$$

References[LINK]

Candanedo, L.M., Athienitis, A., Park, K.W., 2011. Convective heat transfer coefficients in a building-integrated photovoltaic/thermal system. Journal of Solar Energy Engineering. March 2011, Vol. 133

Delisle, V., Kummert, M., 2014. A novel approach to compare building-integrated photovoltaics/thermal air collectors to side-by-side PV modules and solar thermal collectors. Solar Energy. February 2014, Vol. 100, p. 50.

Duffie, J.A., Beckman, W.A., 2006. Solar engineering of thermal processes. 3rd edition, John Willey & Sons Inc., Hoboken, NJ

Gorman, J.M., Sparrow, E.M., Katz, S.D.M., Minkowycz, W.J., 2019 Convective heat transfer coefficients on all external surfaces of a generic residential building in crossflow. Numerical Heat Transfer, Part A:Applications, 2019, Vol. 75, p. 71

Holland, K.G.T., Unny, T.E., Raithby, G.D., Konicek, L., 1976 Free convective heat transfer across inclined air layers. Journal of Heat Transfer, May 1976, Vol 98, p. 189

Incropera, F.P., Dewitt, D.P., Bergman, T.L., Lavine, A.S., 2007. Fundamentals of heat and mass transfer. 6th edition, John Willey & Sons Inc., Hoboken, NJ

Sparrow, E.M., Nelson, J.S., Tao, W.Q., 1982 Effect of leeward orientation, adiabatic framing surfaces, and eaves on solar-collector-related heat transfer coefficients. Solar Energy, 2019, Vol. 29, p. 33