Mathematical model of flat plate solar thermal collector and its validation

. In order to predict solar thermal collector’s performance and optimization of control algorithms, proper mathematical models are necessary. Computer calculation technique provides tools for determination of modern materials impact on improvement of heat transfer inside the collector and minimization of heat loss. Such analysis is impossible by using standard technical datasheet provided by producer or by using empirical formulas. In the paper the authors present a mathematical model of a flat plate solar thermal collector based on the Hottel-Whiller-Bliss equation and criterial formula. The iterative algorithm solved steady state heat transfer equations for a glazed and an unglazed collector. The validation experiment was conducted under class AAA sun simulator for different inlet temperature and solar irradiation values. The unglazed PV/T and glazed solar thermal collectors were tested. For PV/T the relative difference, between measured and computed outlet temperatures, was below 5% and the highest value was reached for the lowest inlet temperature. The validation study showed that the experimental results reached good agreement with simulation predictions. Presented computation algorithm enables to predict influence of geometry changes on collector performance. The model could be used for optimization of the construction without using CFD methods, which need large computation resources.


Abstract.
In order to predict solar thermal collector's performance and optimization of control algorithms, proper mathematical models are necessary. Computer calculation technique provides tools for determination of modern materials impact on improvement of heat transfer inside the collector and minimization of heat loss. Such analysis is impossible by using standard technical datasheet provided by producer or by using empirical formulas. In the paper the authors present a mathematical model of a flat plate solar thermal collector based on the Hottel-Whiller-Bliss equation and criterial formula. The iterative algorithm solved steady state heat transfer equations for a glazed and an unglazed collector. The validation experiment was conducted under class AAA sun simulator for different inlet temperature and solar irradiation values. The unglazed PV/T and glazed solar thermal collectors were tested. For PV/T the relative difference, between measured and computed outlet temperatures, was below 5% and the highest value was reached for the lowest inlet temperature. The validation study showed that the experimental results reached good agreement with simulation predictions. Presented computation algorithm enables to predict influence of geometry changes on collector performance. The model could be used for optimization of the construction without using CFD methods, which need large computation resources.

Introduction
Flat plate collectors are used widely for domestic hot water preparation and domestic heating. Mathematical models of solar thermal collector are used for simulation or optimization of system. Many researchers used the Hottel-Bliss-Whiller equation to compute collectors performance [1]- [4]. Some work applied the Hottel-Bliss-Whiller equation to multiflow numerical analyses [5], air collectors [6] and concentrated collectors [7] have been published so far. The process of the heat transfer between the absorber and an ambient is complicated, especially for the top losses. Heat is transferred by a convection and by a radiation. Empirical formula for the heat loss coefficient was proposed by Klein and further it was developed by others [4], [8]- [10]. The Klein formula has the following form [8], [9]: (1) The Klein formula results in large errors due to calculation of convection and radiation parallel with unknown glass temperature [9]. A mathematical model of solar collector based on analogy between heat transfer and electrical circuit was proposed [2]. The model was developed and empirical equation for glass temperature was presented [9]. Many experiments of convection have been conducted and equations for heat transfer coefficient have been proposed [2], [6], [11]- [13]. The formulas compute heat transfer coefficient as linear function of air velocity, for example McAdams which is a popular formula [2], [6]: During last years some research was conducted to optimize collectors area, storage tank volume and system pattern [4], [14]- [16]. Mathematical model of solar thermal collector used for optimization based on characteristic equation was presented in literature [14], [15]. The Hottel-Bliss-Whiller equation and empirical formula was also used in optimization for calculations of top heat loss coefficient [4]. Optimization of absorber plate geometry was done in reference [17], while complex optimization of geometry, physical properties of materials and system were not reported in the literature. Proposed mathematical model gives opportunity to analyse influence of material properties and geometry on collector efficiency. Moreover, it could be applied in order to reduced collector mass and production cost by using new materials.

Mathematical model
The mathematical model, presented by the authors, based on the Hottel-Whiller-Bliss equation has the following form [1]- [3]: Heat loss coefficient is a sum of bottom heat loss and top heat loss coefficients. The convection from the back of the collector was neglected. The back heat loss coefficient was computed from thermal resistance of the insulation. The top heat loss coefficient was computed from connection of thermal resistances of radiation and convection according to equation [2], [3], [9]: Radiative heat transfer coefficient between plate and glass cover was computed by [2], [3], [9]: For radiation from glass surface to the sky [2]: h r,g-sky = σ ε g (T sky 2 + T g 2 )(T sky + T g )/(T g -T a ) (9) The convective heat transfer in the air gap was calculated from the Hollands equation [2], [3], [ From the glass surface heat is transferred to the ambient by a free or forced convection. For a free convection the Nusselt number was calculated by [2]: The Nusselt number for a forced convection was computed by [2]: Nu = 0.86 Re 1/2 Pr 1/3 (12) The heat removal factor was calculated with a procedure applied by many researchers [2], [3], [18], [19]: where F is fin efficiency.
The heat transfer coefficient in plate channels was computed from the Michiejew equation [20]:  (15) For a computation of the collector's outlet temperature inlet parameters are necessary. The model needs the meteorological and the flow data, such as: irradiance, wind speed, ambient and sky temperature, mass flow through the collector and inlet temperature.

Sensitivity analysis
The influence of the input parameters to the model was analysed by simulations. Testing data sets were prepared for each input parameter (solar irradiation, inlet temperature, ambient temperature, wind speed, mass flow rate). The outlet temperature was calculated for variation of one parameter when others were constant. Procedure was repeated with changed constant parameter. The influence was calculated from equation: where X is the examined parameter, j is the index of data. The results of analysis were presented in the figures below, except for the wind speed, which had the least influence. Influence of inlet temperature locally rising and then falling down. The rise is connected with the working point where efficiency of the collector is close to the zero value. After changing to the negative value influence of the inlet temperature is falling down (Fig. 2).

Convergence and iteration influence
The calculation algorithm had two convergence conditions: • difference between new and old outlet temperature value lower than 0.001°K, • conservation of the energy with an accuracy of 70 J.
The programme reached required convergence in about 10 iteration during sensitivity analysis. Results of convergence and number of iteration were presented in Fig. 8 and Fig. 9 for datasets used for inlet temperature and irradiation sensitivity analysis.
Number of iteration had approximately exponential influence on the convergence (Fig. 10).

Validation experiment
The experiment was conducted under class AAA solar simulator -LA150200 (Eternal Sun) in controllable conditions. The following parameters were measured during the experiment: • irradiation, by SR11 Hukselflux pyrometer, • infrared irradiation, by IR02 Hukselflux pyrgeometer, • temperature, by PT100 A class resistant thermometers, • mass flow, by Optimass 6400C flowmeter.
Weber Sol ECO 2,0 and Ensol 2,0 PV/T solar collectors were examined. Figure 11 presents deviation between measured and computed values of outlet temperature of the PV/T. The experiment was conducted for two different solar irradiations (820 and 1020 W/m 2 ) and for four inlet temperatures (14, 23, 32, 42˚C). Mass flow rate during experiment was set about 1.2 l/min, Temperature in the experimental lab was set between 21 to 29˚C. During each measurement ambient temperature was constant. Air velocity was constant and equal to 1 m/s. Experiment was conducted in steady state, the experiment's data were measured for at least 5 minutes. The uncertainty of the measured data was estimated based on standard deviation and meter accuracy. Temperature's uncertainty was estimated as 0.3˚C, irradiation uncertainty was estimated as 16 W/m 2 . Because of little influence on the computation results uncertainties of wind speed and mass flow rate were not taken into consideration.

Conclusions
The mathematical model of solar thermal collector was presented. The sensitivity analysis showed that inlet temperature value had significant influence on the computation results. 1% change in inlet temperature resulted in about 1% change in outlet temperature. Inlet temperature influence was presented in Fig. 2 and Fig. 3. Second factor with huge impact in sensitivity analysis was solar irradiation (Fig. 4). The influence of mass flow rate depended on the irradiation, ambient temperature and inlet temperature.
The solving algorithm was reaching required accuracy in about 10 iterations. Increasing precision of the solution did not increase significantly the number of iteration (Fig. 10).
The analysis of the comparison between the experiment and simulation leads to the following conclusions. First of all, the presented model achieved a good agreement with the experiment; a maximum difference about 3% was obtained (Fig. 14).
Uncertainty of the measured parameters had not significant influence on the computation results ( Fig. 16 and Fig. 17).
The Klein formula in every simulation underestimated collector efficiency and outlet temperature (Fig. 12, Fig. 13 and Fig 14). Moreover, every use of the Klein formula gave higher difference between experiments than the model. Computation costs of Klein formula algorithm and the model were similar. Both methods needed about 4 iterations during validation simulation. Simulation with Klein formula in most cases reached better convergence but results were not as closed to the experiment as the model's results. The presented model gave more accurate results than a model based on the Klein formula. According to the experiment for PV/T collector the model can be used in order to achieve 5% accuracy. This work was supported by the own research fund of The Szewalski Institute of Fluid-Flow Machinery PASci (FBW27).
One of the authors (M. K.) was supported by the Foundation for Polish Science (FNP).  Fig. 19. Comparasion between models in respective to energy convergence.