Energy Convers.Mgmt Vol. 31, No. 5, pp. 409-418, 1991 Printed in Great Britain. All rights reserved
0196-8904/91 $3.00+ 0.00 Copyright © 1991 PergamonPress plc
THE SOLAR WATER HEATING SYSTEM WITH NATURAL CIRCULATION ASSISTED BY AN AUXILIARY ELECTRIC HEATER--PERFORMANCE MODELING LIN WENXIAN, LU ENRONG and WANG DONGCHENG Solar Energy Research Institute, Yunnan Teachers University, Kunming, Yunnan Province, People's Republic of China 650092 (Received 7 October 1988; received for publication 20 April 1990) Almtraet--This paper is the first part of the research on the performance of a solar water heating system with natural circulation assisted by an auxiliary electric heater. Based on the up-to-date best models that simulate the performance of such a system--Morrison-Braun's fixed-node nonconvection-mixingmodel and variable volume fully-convected-mixing model, this paper combines the advantages of these two models and abandons their disadvantages and develops a more precise and widely-used fixed-node equally convected-mixing model and also includes the effects of the auxiliary electric heater. The experimental results verify the developed model's validity and its precision. Solar water heating system modeling Temperature convected-mixing model
Natural circulation Auxiliary electric heater Performance distribution Finite-difference method Fixed-node equally Solar energy saving ratio
NOMENCLATURE
Collector aperture area (m2) Storage tank side area (m2) A x = Storage tank bottom area (m2) bo, bou), Bo, ~ , B. = Decimal part of numbers Specific heat (J/kg °C) d = Connecting pipe inner diameter (m) d,= Collector risen pipe inner diameter (m) E = Node of electric heater in storage tank Electric heater power (W) f = Friction coefficient (dimensionless) f,= System solar energy saving ratio (dimensionless) F = Control factor (dimensionless) F ' = Collector efficiency factor (dimensionless) Collector heat removal factor (dimensionless) g i = Control factor (dimensionless) h T . - b = Thermosyphon head between a and b in a pipe (N/m) H , = Storage tank height (m) It= Solar radiation intensity ( J / m 2 s) 1 H = Node of inlet hot water coming from collectors into storage tank L i = Pipe length (m) System mass flow rate (g/s) K~I = Mass flow rate removed by load (g/s) Nc,N,,N. ffi Integral part of numbers Re = Reynolds number (dimensionless) r.= Atmosphere temperature (°C) rR= Control temperature (°C) ~'L= Collector total heat loss coefficient (W/m2 °C) Pipe heat loss coefficient (W/m2 °C) Storage tank heat loss coefficient (W/m2 °C) A N = Height difference between tank and collector (m) At ffi Time interval (s) ~ x - - Node interval (m) Collector slope (deg) O= Solar incidence angle (deg) p = Density (kg/m3) Kinetic viscosity (kg/m s) Bond's friction coefficient (dimensionless) rid = Collector average daily efficiency (dimensionless) Collector transient efficiency (dimensionless) Ac = As=
E C M 31/5--A
409
410
WENXIAN et al.:
SOLAR WATER HEATING SYSTEM
INTRODUCTION It is very difficult to theoretically analyse the performance of a solar water heating system with natural circulation because of the system's own intrinsic operational characteristics. Although many theoretical analysis models having different complexities and precisions [1, 2] have been put out, the best models, at present, are the ones developed by Morrison, Braun and others [3-9]. They simplified the complexity of the theoretical analysis of the system performance with reasonable treating techniques and obtained the up-to-date best models of the fixed-node nonconvectionmixing model and the variable volume fully-convected-mixing model. These two models' accuracy is not so high because they correspond, respectively, to the large and small mass flow rates, and their application limits are narrow. Therefore, it is necessary to improve these two models and to develop a more precise, and widely-used theoretical analysis model. In a solar water heating system, the auxiliary heating is a very important part, but there is little work that has been done. Heretofore, there are no detailed research reports of the effects of the auxiliary heater on the solar water heating system performance. This is because there is still no ideal analysis model, so far, and the auxiliary heating increases the difficulty of analysis. The model developed in the present work includes the part of the auxiliary heating, so we can make the detailed analysis and research on the effects of the auxiliary electric heaters. THEORETICAL ANALYSIS MODEL The performance simulation is carried out with the finite-difference method. System Temperature Distribution Water temperature distribution in the collector
Along the water flow direction in the collector, the collector can be divided into N~ equal parts (nodes). The water temperature distribution is represented by the node temperature distribution. The temperature Tk of the node K (1 ~
(1)
where
Ta, = Ta + 4F,,(~=).IFRUL, F,
=
F'ULA~/(mcCpN=)
and F'UL can be calculated by: F ' U L = _ mtCp In (1 Ac
\
FR ULAe.'~
mtC. ]"
(2)
So, the total useful energy gained by the collector is: Qu = rAcAt [FR(T~t).It - FR UL( Tf, i - Ta)]
(3)
where the modification coefficient r is calculated by: r =
FR(rhc) the[1- exp(1 - F'ULAc/rhcCv) ] = V~(,n,) rh,[1 - exp(1 - F'U~Ao/rh, CO]"
(4)
Water temperature distribution in the connecting pipes
As the thermal loss coefficient is generally very small, and the surface that loses heat energy is also very small, we can treat the upper and lower connecting pipes as two nodes. Each node temperature is represented by its respective connecting pipe average water temperature, viz. T,, --- Ta + (/'4
-
Tp2 = Ta + (T2 -
,+oc, [ 1 ,*oc. [ 1 Ta) ~ V,) ~
-
exp exp
(VA),,-] ~-] (VA),:] ~ ]
(5) (6)
WENXIAN et al.: SOLARWATER HEATING SYSTEM
411
and the collector inlet temperature T~ and the hot water inlet temperature T3 on the storage tank are: 71 -- Tf.i = Ta -I- (74 - Ya) exp[-(UA)p~/YncCp]
(7)
73 = T, + (T2 - T , ) e x p [ - ( U A ) p J m o C p ] .
(8)
Water temperature distribution in the storage tank We can divide the water body in the storage tank into N, equal parts (nodes) along the vertical direction. Suppose the auxiliary electric heater is located at the node E, the hot water coming from the collector entered the storage tank at the node IH, the water in the storage tank was out at the first node and entered the collector loop, the load supplied cold water to the storage, tank at the first node (controlled by a floating-sphere valve) and removed the hot water from the storage tank at the node N, (the top of the storage tank). Generally speaking, the node energy differential equations in the storage tank include the nonsteady term (time term), convecting term, thermal diffusion term and thermal source term. If considering simultaneously these four terms, it is very complex and very time-consuming to simulate the system performance. In practice, it is almost impossible to theoretically simulate the system performance if some special treating techniques are not being used [10-12]. The solutions of the energy differential equations in the storage tank can be obtained in two stages: the first stage, not considering the convecting term and solving the equation including the remaining three terms. Such equations can be easily solved. The second stage considers the convecting term effects by using analysis treating methods and finally obtains the water temperature distribution in the storage tank. Temperature distribution T~t(i) not considering the convecting term. In this case, the energy differential equation of node i (1 ~< i ~
(pCp) a~._~(ti)- O-x d I K, ~ ]
- U'[Ts,(i) - T~,] + g,FEw d-~"
(9)
The boundary and initial conditions are:
aTsl Ot
x=O,t = 0 ,
~aT~l ~=",,' = 0
(lO)
T,,(x, O) = TLo.
(1 l)
. The finite-difference general solution of the above equation is (see Appendix)[13]:
a~°) Tsl(i)
= ~-i- - " Ct°-~) osl) . ~,
D/_a-a~O)i+lT~°)risl~+
1) + b
(12)
where
2/1 1 ~ -I al°_) , = ~ x / -\'lxi ~ + ~ Ji - 1] ~i+1 ~"~-'XX
)
/~i+1
a<0)- : , , r ~ Ax
(13)
b = (a~°) - al°)_, - - ,,to) U'x) ~°)(i) + giFEw + U'T, Ax ~i+ I ~ T and
a0 -- 0, aN, + 1= 0. In the above equations, the superscript (0) represents the time interval t - At --, t (not including time t) and is the last time interval; T,~(i) represents the value in the current time interval t ~ t + At
.412
WENXIANet al.: SOLARWATERHEATINGSYSTEM
(not including time t + At). According to the analysis, under normal conditions, At can be large, up to 15 min and has little influence on the system performance if U' is not so large, and this requirement is generally met. Because the specific heat Cp of water has little change with temperature, we can take Cp = 4180J/kg °C in the above equations. If not involving the thermosyphonic term, the water density can be taken as p = 990 kg/m 3, otherwise, we must use the following correlation [14]: p =1000.5-0.0762439T-0.00349823T 2
(14)
K = 0.551 + 0.0026T - 0.0000132T 2.
(15)
and K is calculated by
Temperature distribution Ts(i) considering the convecting term. Convection mixing is caused mainly by the supplement of hot water coming from the collector (refer to the collector flow), the removing of hot water from the storage tank by the load (refer to the load flow) and the fast speed heating of the auxiliary electric heater. When analysing the effects of the convection mixing, we can consider these factors separately and then carry out the algebra treating methods [7]. (a) Temperature distribution Ta(i) after the collector flowing. In this case, the load flow is not considered and only hot water coming from the collector enters the storage tank, and at the tank's bottom, the colder water coming from the tank enters the collector. In the time interval t ---, t + At, the volume of the hot water coming from the collector is V, = Srh¢dt/= moat/. On the treatment of the effects of the convection mixing caused by hot water coming from the collector, Morrison and Braun took two extreme cases: the fixed-node nonconvection-mixing model supposed that the hot water coming from the collector does not convectively mix with the water layer in the storage tank that has the same or the closest temperature as the entering hot water, and it is approximately valid when the mass flow rate is very small. The variable volume fully-convected-mixing model supposed that fully-convected mixing occurs between the hot water coming from the collector and the water layer in the storage tank that has the same or the closest temperature as the entering hot water, and it is approximately valid when the mass flow rate is large enough. Therefore, the accuracy of these two models is not very high because of their supposed condition limits. On the other hand, they are inconvenient to widespread use in practice because of the narrow operating limits. Based on these two models, this paper supposes that the temperature T~(= 7"3) of the hot water coming from the collector is nearest to the temperature T,I(G) of the node G in the storage tank, then the hot water, with volume Vo and temperature To, coming from the collector mixes equally-convectively with the various nodes between nodes 1H and G. It is apparent that such treatment combines, to a great extent, the advantages and abandons the disadvantages of Morrison-Braun's two models, and it is predicted that the model developed in this paper (refer to the fixed-node equally-convected-mixing model) has a higher precision and wider operating limits. It can be used not only at the low mass flow rate but in the high mass flow rate system and especially at the middle mass flow rate (this is the general case). So, it can be a widely used simple model. Temperature distribution T,2(i) after considering the collector flow is: if 1 ~
j--I
(16)
if (Z - Arc- 1) ~< i ~<(Y - n~ - 2), then
f
Ts2[Y-j~|(nc(:)+l)-I 1.....Ts2[Y)j-I ~(nc(j)-.l-l)l=T~2(Y-k k+l Tv.
]
~ (ncf:)+ 1) =bc(k+l)T~(Y-k)+[1 -b~(,+l)]r'a(Y-k - 1) -1
07)
WENXIAN
SOLAR WATER HEATING SYSTEM
et al.:
413
if ( Y - n~ - 1) ~< i ~< Y, then
Ta ..... Te(Y- n,:)= T~(Y) Ta(Y - n¢ -
boT',2(Y)+
I) =
(I -
b¢)T~2(Y-
I)
(18)
if (Y + 1) ~
T~2(i )
(19)
( Tsl(i), otherwise
when Y/> E,
f T,~z,F= l,i = I ,y f +.,Ns
T,2(i) = "~ Ts1(i), r = O, i -
,Ns
(20)
[ Ta(i),F=O,i=E,...', Y."
(b) Temperature distribution T~(i) after considering the load flow. In this case, only hot water is removed from the storage tank top and cold water is supplied at the bottom during the interval t. The removed hot water volume by the load is Vi = m~At/p. Set N1 + Bl = Vl/Ux = mlAt/(pAxAx), then all node temperature distributions T~(i) (1 ~
Ts(i)= BiTs2(i- Nl -- I) + if V~> Vv, set N n + B n = N I + B I - N c - B c ,
{
T,(l)
=
T,(2)
....
(I -- BO
Ts2(i- Nt)
(21)
then
T,(N.)= T,.o
=
T,(N.+ 1)= ~B'TL°+BcT''(1)+(1 - B n - Bc)T~2(I), (B~< 1 - B n ) (B, TLo+ (1 -- B,) T~,(1), (Be/> 1 - B,)
(22)
for N, + 2 <<.i <<.N,, T,(i) is calculated with equation (21). It can be obtained that the water temperature T4 coming from the tank bottom and entering the collector loop is:
T4 = (N,+ B,)TLo+ ~, Ts(i)+ BcT~(N¢+ I)
(N~+ B,+ N~+ Bc)
(23)
i=l
and the outlet water temperature To coming from the tank top is:
To= I li
~
Ts(i)+BIT,(Ns-N,)I/(N,+B,).
(24)
ffi N s - N I + 1
So far, the temperature distribution of the system loop can be obtained.
Calculation of the System Mass Flow Rate System thermosyphonic head Ahr AhT is calculated by [1, 4] Ahr = hTp -
hTm
(25)
and I
hTp =
In / H i~l p,~Axg -- Axg(ps~ + P~m) 2 + P4-~g( ~+ AH + Ax/2)
~hTm = HIgj~ PcflN~ + p2- 3g[AH + (IH - ½)Ax] .
(26)
WENXIAN et aL: SOLAR WATER HEATING SYSTEM
414 System flow friction AP
AP = APy + APj
(27)
where APy is the friction caused by the water circulating in a loop through the collector, the connecting pipes and the storage tank. For the solar water heating system with natural circulation, the water flowing in the pipe is generally in the laminar flow (Re ~< 2300). So, we can c a l c u l a t e f w i t h f = 64/Re (if Re > 2300, take Re = 2300 [4-7]). Considering the non-fully-developed laminar region in the pipe, there is a friction coefficient modification factor M in f [6, 7]: f = M64/Re
(28)
L ~-0.96 M = 1.0 + 0.038 \D---?-.~Rej
(29)
where
Hence, AP e is calculated by 128rh¢ FIt, L_2LM, It2L2M2 It3L3M3 (/t4M4 + bqMs)L4] L d4 + d---4--+ nd---~-l -d~ J"
(30)
.~" 2 S M 8en~ ~ p d 4 ~ ' ~,
(31)
a P ~ = - - rcp
APj is calculated by
aP;
=
The values of ~i of the various elements can be found in the relevant books. System mass flow rate rhc To maintain the system's natural operation, AhT must be no less than AP. Under the condition of AhT = AP, we can obtain the equation to calculate rh0: rh¢ = [x/(BC) 2 + 4AQ - BC]/2/A
(32)
where /1 - n ,.~ pidt4 B = 128/d4/n
C = It, L~Mi + It2L2M2 + l~3L3M3(d/d~ )4In + (it4M4 + #sMs)L4 Q=g
Axe'p,,-
pc; N~ - A x g ( p , , + p a n ) / 2
i-1
"=
+ P4-~g (HI + A H + AX/2) - P2-3g [AH + (IH - ½)Ax] and/z is calculated by [14] It = 0.164323 x 10 -2 - 0.393398 x 10-4T + 0.43606 x 10-6T 2 - 0.180044 x 10-ST 3.
(33)
System Performance Parameters (a) Collector transient efficiency ~/t t/, = rh¢Cp(To - T~)I(ITA¢).
(34)
(b) Collector average daily efficiency t/a [9] r/d = rhcCp(To - T~) dtllrA, dt = Y.mcCp(To - Ti)IZITA¢.
(35)
(c) System solar energy saving ratio f~ [5-7] f~ = 1 -- QJ(Qt4 + QL,)
(36)
WEN-XIAN et al.:
SOLAR WATER HEATING SYSTEM
415
0.3 -
0.2 E 4=' O
3= o.1
g
"I-
I
I
I
1
9:00
11:00
13:00
15:00
17:00
Time
Fig. l. The distribution curve of hot water removed (%) by the load.
where Q,, QLd and QLs are the consumed auxiliary electric energy, the useful energy removed by the load and the thermal energy lost by the storage tank in the duration of operation in a day, respectively. Solution Procedure
First, the system's temperature distribution is calculated with the mass flow rate in the last time interval and data of solar radiation and atmosphere temperature in the present time interval. Secondly, the mass flow rate in the present time interval is calculated. Thirdly, this procedure is 1,00
"°° I
..... ----+
Model. 1
j
ModeL 2 Model 3
..1..._~
I
Measured
'"
ModeL 1
.....
Model 2
----t-
Model 3 Measured
rJ
0.7,5
0.75
LX,' i
rj
1- rJ I
-r
J
0..'50 14:00
"r
12:00
rJ
i 0'25i
0.25
o 20
30
40
50
60
T (°C)
Fig. 2. The transient temperature distribution in the tank (12:00).
20
iL I
30
40
I
I
50
60
T(°C)
Fig. 3. The transient temperature distribution in the tank (14:00).
WENXIAN et al.: SOLAR WATER HEATING SYSTEM
416 90
CaLcuLated with model 1
80
Measured
+
70
CoLLector outlet
6O ~o
40' 30
CoLLector inlet
~_
=I*
20'
~.
~"
+
+
=I"
1.
t"
"r
lOo lO:O0
.I
I 12:00
11:00
I 13:00
I 14:00 Time
I 15:00
I 16:00
I 17:00
Fig. 4. The collector inlet and outlet temperatures.
repeated with the mass flow rate in the present time interval and the data of solar radiation and atmosphere temperature in the next time interval. EXPERIMENTAL
VERIFICATION
Experimental Device Some system parameters of the experimental device are: Ac= 1.872m 2, FR(Z=),=0.794, F~UL = 8.38 W / m 2 °C, K~= = 1.004 - 0.182(1/cos0 - 1), mt = 0.02 kg/s, dl = 0.015 m, d = 0.032, n = 9. The connecting pipes have an inner diameter of d = 0.032 m, and the pipe thickness (black rubber pipe) is 0.005m. K = 0 . 1 6 W / m ° C , L~ =2.65, L2= 1.14m. The storage tank has the dimensions 0.7 x 0.35 x 1.0 m 3. Ns = 20, E , = 1 kW, E = 18, IH = 16, x = 0.0 m, TR = 55°C, H = 0.4m. There are 10 pairs of copper-brass thermocouples in the storage tank. The mass flow rate is measured with the ink tracing method. For the sake o f simplicity, the simple curve shown in Fig. 1 was used in the experiment. The load removed 200 kg of hot water from the storage tank between 9: 00 and 17: 00.
Experimental and Theoretical Simulating Results The experimental and theoretical simulating results are shown in Figs 2-6. To judge conveniently, the theoretical simulating results with the two Morrison-Braun models (models 2 and 3) are also shown (model 1 is the one developed in this paper). 16
14 12 10 8
!Y
_
6. 4
+
Measured
20 10:00
I
I
I
I
I
I
I
11:00
12:00
13:00
14:00
15:00
16:00
17:00
Time
Fig. 5.
The
system thermosyphonic mass flow rate.
WENXIAN et aL: SOLAR WATER HEATING SYSTEM
417
80-
6O 50 40 30"
÷
CaLcuLated with model I
20 --
+
Measured
10-0 10:00
I
I
I
I
I
I
I
11:00
12:00
13:00
14:00
15:00
16:00
17:00
Time
Fig. 6. The collector transient efficiency.
From these figures, it can be seen that the theoretical simulating results of model 1 agree better with the experimental results than the two other models. The experimental results verified the developed model's validity and its precision.
CONCLUSION
The model developed in this paper (fixed-node equally convected-mixing model) is the best model at present; it is widely applicable and it can be used to simulate the performances of various complex solar water heating systems with natural circulation.
REFERENCES 1. M. S. Sodha, S. S. Mathur and M. A. S. Malik (Eds), Reviews of Renewable Energy Resources, Vol. 1. Wiley, New York (1983). 2. W. C. Dickinson and P. N. Cheremisinoff (Eds), Solar Energy Technology Handbook (Part B). Marcel Dekker, New York (1980). 3. G. L. Morrison and D. B. J. Ranatunga, Sol. Energy 24, 55 (1980). 4. G. L. Morrison and D. B. J. Ranatunga, Sol. Energy 24, 191 (1980). 5. G. L. Morrison and C. M. Sapsford, Sol. Energy 30, 341 (1983). 6. G. L. Morrison and H. N. Tran, Sol. Energy 33, 515 (1984). 7. G. L. Morrison and J. E. Braun, Sol. Energy 34, 389 (1985). 8. G. L. Morrison, Sol. Energy 36, 377 (1986). 9. R. Uhlemann and N. K. Bansal, Sol. Energy 34, 317 (1985). 10. A. Cabelli, Sol. Energy 19, 45 (1977). 11. K. DenBraven, ASME J. Sol. Energy Engng 108, 105 (1986). 12. A. M. C. Chan, P. S. Smereka and D. Giusti, ASME J. Sol. Energy Engng 105, 246 (1983). 13. S. V. Patankar, Numerical Heat Transfer and Fluid Flow. McGraw-Hill, New York (1980). 14. M. Vaxman and M. Sokolov, Sol. Energy 37, 323 (1986).
APPENDIX In the case of not considering the convecting term, the energy differential equation and the boundary and initial conditions are equations (9)-(11), where U ' = (UA)s/(H~Ax), and F = { O , ifrhl=O
or
TN(Ns)TR
I, otherwise g~ ~
0, i ~ E 1, i =E.
From Ref. [13], the finite-difference general solution of the above equation is equation (12).
W E N X I A N et al.:
418
SOLAR W A T E R H E A T I N G SYSTEM
The solutions of node 1 and node N, are determined with the boundary conditions: (a) Node 1 The energy differential equation for node 1 is:
Coc~),dT,,(1) at
d
=~S
g~ OT, i(1)] -- U'[T,t(1 ) - Ta] + -g,FEw -
~j
A,x
Integrating in the control volume of node 1, we can obtain: ~xx ( - q ) d t dx - U'ataxT{l)(l)+pg, &x + u'T, Atax
(pCp)I[T'~(1)- ~°)(1)]--Jb Jt as qb
:
0, qm : ~
+
[T~I)(1) -- ~0)(2)1 -----a~°)[~l)(1) -- ~0)(2)1
hence.
al°iT, i ( 1) = a~°) ~°i)(2) + [at °) _ at °) _ V'~:] T~°?(I ) + FII,FE,, +
L A,x
u'r,/ax.-I J
(b) Node N s The energy differential equation for node N, is: (~" OTsi(Ns) ~ r
P%~"
~
dT, I(N,)- ]
-- o~ L X . , ~ j -
g~FEw
u'[r,,(N,)- r,] + A,--S-
Integrating in the control volume of node N,. we can obtain:
r.r '-'o~x (-q)
(pCp)~,[T,~(N,) - T~°)(N,)iAx : ,J- ,/t as qu =
O,
qm = - ~
+
dt dx - U'AxAtT~°)(N,) + L - ~ , x [T~°)(N, -
I) - T(,°)(N,)I
hence.
a(°)T~, ,l,~-~ ~ - a(°)-N, lTl°)(N,-1)
[a~,°)_ a ~I°)_ 1 - U'AxIT
Combining the above equations, we can obtain equations (12)-(13).