Distributed Parameters Dynamic Model of a Solar Fresnel Collector Field

Distributed Parameters Dynamic Model of aSolar Fresnel Collector Field ⋆ Marco Spoladore ∗Eduardo F. Camacho ∗M. Elena Valcher ∗∗ ∗System and Automation Engineering Department, University of Seville, Spain (e-mail: [email protected],eduar[email protected]) ∗∗ Department of Information Engineering, University of Padova, Italy (e-mail: [email protected]nipd.it) Abstract: This paper describes adynamic model for alinear Fresnel collector field of asolar refrigeration plant. The collector field concentrates solar radiations on atubethat heats up water that is used byan absorption machine to produce chilled air for refrigeration purposes. The model takes into accountthe solar radiation losses produced bythe mirrors and the absorbing tubestructure as well as the temperatures in 64 segments of the receiving tubeand water. Although the dynamic model is of high dimension andnonlinear, the identification problem can beexpressed as alinear problem in model parameters. The least square methodwas used for identification. The model was validated comparing the simulation results with the plantdata for differentoperatingconditions. 1. INTRODUCTION In today’s climate of growing energy needs and increasing environmental concerns, alternatives to the use of nonrenewable and polluting fossil fuels need to besearched for. One suchalternativeis Solar Energy.The sun creates its energy through athermonuclear processthat creates heat and electromagnetic radiation. The heat remains in the sun and is instrumental in maintaining the thermonuclear reaction. The electromagnetic radiation streams out into space in all directions. Only avery small fraction of the total radiation reaches the Earth, where it is collected. Solar energy is radiated from the sun through forms of heat and light, and all Solar Technology Systems capture the energy of the sun byabsorbing lightas heat. There are two main types of solar technologies: Thermosolar technology and Photovoltaic technology.Thermosolar Power Systems concentrate sunlight, usually with mirrors, to heat to high temperatures afluid that drives an engine. This approachdiffers from the one used byPhotovoltaic Solar Power Systems, where lightinteracts directly with special materials to separate charges and generate electricity.The world of Solar technologies consists of 4major areas: (1) Concentrating Solar Power Systems; (2) Photovoltaic Systems; (3) Solar Heating Systemsand (4) Solar Lighting Systems. This paper focus on Concentrating Solar Power Systems, in particular linear concentrator systems in whichreflectivematerials focus the sun’s heat energy that then drives agenerator to produce electricity.Specifically,wewill consider the Linear Fresnel Reflector System of the Seville University’s Thermosolar Plant, that will bedescribed in the next section. The paper provides adynamicmodel for this system whose main purpose is to help developing a ⋆This work was supported in part bythe Spanish Ministry of Science and Innovation,Junta de Andalucia and Gas Natural SDG S.A. under grants DPI2008-05818, P07-TEP-02720 and PI-0709/2009 control scheme for this typeof plants.Fresnelconcentrators based on the approximation of the optical behavior byray-tracing havebeen developed (Larbi et al. [2000], Facao et al. [2009]). These non-dynamical models were developed to evaluate and optimize collectors design. This paper presents adynamic model of the Fresnel collector field that allows to simulate optical factors suchas the shades produced bythe facets, the tubestructure, the temperature of the metal tubeand of the heated water. Adynamic model of parabolic trough collectors has been developed (Camacho et al. [1997]) for the solar platform of Almeria (PSA). The paper is organized as follows: section 2describesthe Linear Fresnel Collector System and section 3provides the related optical model. The heat interchange in the receptor pipeis described bythe thermal model in section 4, while model validation results are shown in section 5. 2. LINEAR FRESNEL REFLECTOR SYSTEM Linear Concentrating Solar Power(CSP) collectors are one of the three types of CSP systems in use today.They include the twomajor types of linear concentrator systems: Parabolic Trough Systems and Linear Fresnel Reflector Systems.Linear CSP collectors capture the sun’s energy with large mirrors that reflect and focus the sunlight onto alinear receiver tube. The receiver contains afluid that is heated bythe sunlight, and then used to create superheated steam that spins aturbine that drives a generator to produce electricity.Alternatively,steam can begenerated directly in the solar field, eliminating the need for costly heat exchangers. Linear concentrating collector fields consist of alarge number of collectors in parallel rows that are typically aligned in anorth-to-south orientation to maximize both annual and summertime energy collection. With asingleaxis sun-tracking system,this configuration enables the Proceedings of the 18th World Congress The International Federation of Automatic Control Milano (Italy) August 28 - September 2, 2011 978-3-902661-93-7/11/$20.00 © 2011 IFAC 14784 10.3182/20110828-6-IT-1002.02992 mirrors to trackthe sun from east to west during the day,ensuring that the sun reflects continuously onto the receiver tubes. In linear Fresnel reflector systems, flat or slightly curved mirrors mounted on trackers on theground are configured to reflect sunlightonto areceiver tubelocated above these mirrors. Asmall parabolic mirror is sometimes added atop the receiver, to further focus the sunlight. In the thermosolar plantsituated in Seville’s Engineering College, the solar energy collected byaLinear Fresnel collector field is used byan absorption machine to provide cooling for the College. Fig. 1. Seville’s Thermosolar Plant The solar plantfunctioning can beschematically divided into three phases: (1) Collector mirrors reflect sunlightto an absorber pipe situated on aplane parallel to that of collector, but located over it. (2) Aheat-transport fluid (water) flows through the pipe, absorbing heat from the walls of the duct to warm up, increasing its own temperature. (3) The hot fluid enters the absorption machine, that exploits water’s heat to refrigerate the buildings of Seville Engineering’s College. The thermosolar part of the plantconsists of: •Fresnel flat collectors. Infact, they are not completely flat, but slightly curved,with acurvature obtained by aprocess of cold bending. •Areceiver, formed byasingle pipe, with covering. The receiver is situatedon aplane parallel to flat collectors, over them. Inside the pipeflows the heat transport (or heat-carrier) fluid. •Areceiver carrying structure. •Asolar tracking mechanism. In this concentrator typology,eachmirrors rowmoves independently. The main advantages of acollector field of this type are: mirrors and tracking system’s lowcost; almost flat collectors placed next to the ground, so that wind loads are reduced; efficientground use because collector rows can be placed side byside; easy access to mechanic moving parts and surfaces. In Table 1technical data and the principal rating values of the elements that constitute the solar plantare reported. Occupied Ground Extension 480 m2 Total Reflecting Surface 352 m2 PlantOrientation East-West Nr. of Receptors (Nr. of pipes) 1 Receptor’s Length 64 m Receptor Typology Cavityreceptor with secondary reflector and glass cover Receptor Plane’s Height4mover collectors plane Receptor width 0.3 m Absorber Typology Steel pipeDNI 1.4541(AISI 321): stabilized austenitic stainless steel Heat-Transfer Fluid Water Duct Nominal Pressure 13 bar Reflecting Mirror Rows Number 11, eachone divided into 2middle row Single ReflectiveMirror Length 4m Single ReflectiveMirror Width 0.5 m Total Number Of ReflectiveMirrors 176 Mirror Nominal Reflectivity0.92 Optical Concentration Ratio25 Table 1. Technical data of the solar plant. 3. OPTICAL MODEL The optical model developed has general validityand applies to anysolar plantof this typology.In the optical model, reflectivemirrors are considered to beflat, because of their high radius of curvature, and theorigin of the reference system is placed at the beginning of the mirrors’ middle row: the resulting model has the inputs and the outputs indicated in Figure 2. The elaboration of all input data provides all the magnitudes necessary to monitor the plant’s functioning and for the thermal model’s setting. Fig. 2. Optical Model Scheme. In order to understand the magnitudes of the model, a single mirror is represented in Figure 3, with all the parameters necessary to perform calculations on the model. Referring to Figure 3, ori is the orientation angle between the xaxis and the north-south direction, the vector −→ i 18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011 14785 is the solar incidencevector that defines the waythe sunrays affect eachmirror, and its components depend on the values of solar altitude,azimuth and orientation (see Spoladore [2010] for amore detailed description). Fig. 3. Reference system applied to asingle mirror. The optical model is obtained bycombining atwodimensional model and athree-dimensional one. The 2D- Model is developed to calculate the inclination of each mirror row, since eachrowhas the same tilt throughout its length. The model is then used to calculate shadows effects connected with mirrors’ tilt and it is the basis for the 3D-Model. First, the position of eachmirror rowin an ordinary dayat anytime has to becalculated. Eachrowof mirrors has an independentsun-tracking systemthat rotates along the z axis, varying its inclination according to the solar altitude, or elevation.The solar incidence vector is projected on the xy-plane to define its position with respect to the mirrors’ tilt. The projection is again avector, −→ i2D,that still depends on the azimuth value α.The expression for the two-dimensional vector is the following one, where φand θsindicate orientation angle and elevation, respectively: (1) α≥0:−→ i2D=  cos θs·cos (φ+α) sin θs  (2) α<0and |α|≤φ:−→ i2D=  cos θs·cos (φ−|α|) sin θs  (3) α<0and |α|>φ:−→ i2D=  cos θs·cos (|α|−φ) sin θs  The vector −→ i2Drepresents the sunlightthat affects the mirror while −−→ r2Ddenotes the lightreflection on it. Once the incidence vector is determined, the direction of −−→ r2D is forced to bethe same as the receptor’s, to warm it. Clearly the value of the reflection vector is calculated for eachmirror row. Expressing the coordinates of eachrowas (xm,ym), where xm,ymare the rowpositions on xand y axis, the reflection vector is obtained through theformula: −−→ r2D="r2Dx r2Dy#=     XR−xm p(XR−xm)2+(YR−ym)2 YR−ym p(XR−xm)2+(YR−ym)2      where (XR,YR)are the receptor coordinates. Once the components of −−→ r2Dare known, the angle β,formed by−−→ r2D and the horizontal x-axis, can becalculated as: β=arctan r2Dy |r2Dx|. The only twomagnitudes that are still to bedetermined are the angles αand θ:the vector scalar product and laws of reflection can beused to calculate them. Mirrors rowinclination can becomputed, differentiating the calculation based on the various possible cases. xm<0,γ≥β:ρ=π/2−β−θ xm<0,γ<β:ρ=π/2−β+θ xm=0:ρ=θ xm>0:ρ=β+θ−π/2. With this procedure, the inclination value of eachmirrors row, at every hour of every day,can bedetermined. 3.1 Mirrors Shadows The 2D-model is utilized also to calculate the total shadows on thecollectors surface area. Total shadows that affect mirrors area are the sum of twodifferentcomponents: •shadows projected byeachmirror rowon the previous one; •shadowproduced bythe receptor on the total collectors area. For the particular path of the sun abovethe plantthe shadowproduced byeachmirror can only befrom a mirror onto the previous one. The shadowprojected by the receptor on the collector’s surface is usually smaller in percentage than the mirrors’ total shadowand it affects only the negativexsemiaxis. Once the total shadowin the 2D-model has beencalculated, the total shadowsurface is computed bymultiplying the shadowwidth value along the xaxis bythe zaxis length, taking into accountthe orientation of the zaxis with respect to the south and the mirrors currentinclination. The proper mirror shadowis the sum of the shadows that eachmirror projects on the previous one. The possible projected shadowdepends not only on the sun elevation, but aboveall on the corresponding tilt of the adjacent mirror axes. Therefore formulas for proper mirrors shadow calculation are differentiated accordingto mirror positive or negativeinclination. Toperform the calculation trigonometric relations are used. Also the receptor shadow, that is the shadowproduced bythe receptor on the mirror’s reflectivesurface, is calculated using trigonometric relations. The proper shadow proportion pps is theratio between the proper mirrors shaded area and the total collector area (352 m2), while the proper shadow factor fps is the percentage of mirrors’ available reflectivearea, after 18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011 14786 removing proper mirrors shadow. Clearly,fps =1−pps. The receptor shadowproportion,prs,and the receptor shadow factor,frsare similarly defined and in turn are related as frs=1−prs.Then, sum proper mirrors shadowand receptor shadow, the total shaded zone on the reflectivesurface is obtained. ptot =pps +prs 352 ftot =1−ptot. 3.2 Optical Losses The real radiation that reaches the receptor’s pipeis not the same one that affects the reflectivemirrors, because of optical losses. The incidentsolar radiation is fully exploited only if its direction is perpendicular to the reflectivemirror’s surface, but this never happens because mirrors are tilted to reflect solar radiation on the absorber’s pipe. Thus, there is always aloss of direct radiation that depends on the angle κthat asingle mirror formswith the perpendicular to the direction of the solar inclination. The angle κis determined in twodifferentways, depending on the mirror’s inclination value ρ: ρ>0→κ=π/2−γ−ρ ρ<0→κ=π/2+ρ−γ. Optical losses factor of the i-th mirrors rowisdefined as foli=cos κi,where κiis the angle of optical loss relativeto i-th row. Optical losses factors willbeutilized tocalculate the real solar radiation that affects receptor’s pipe. Nowthe 3D-Model is used to calculate the real solar radiation affecting the receptor’s pipe. In order to compute this magnitude, the receptor’s unlit part has to beascertained. The unlit part depends on the azimuthal solar component, that influences solar radiation affecting receptor, therefore a3D-Model has to beutilized. All three components of the solar incidence vector are considered, with the zcomponentvarying between 0and 64 meters, that is the zaxis length. The calculation is performed for eachmirrors row, and hence 11 differentvalues are obtained. So, the value of unlit part δupAis the average of the 11 values: δupA=1 11 11 X i=1 δupi. The receptor unlit part factor fup,that represents a correction factor used in the following calculation of real radiation, is given by fup =64 −δupA 64 . All calculations are made byassuming that the total sunlightreflected bythe mirrors affects the absorber’s pipe or secondary reflector. This is guaranteed bythe large mirror’s radius of curvature. In fact, due to high mirror’s radius of curvature there are no reflected radiation losses. The real solar radiation that affects the receptor is different from the theoretical one, because it takes into accountall the losses and correction factors obtained in the previous subsection. The contribution that eachmirror gives to the solar radiation that affects the receptor (in Watts /m2)is computed according to the following equation: radreali=dr ·Ai·cos κi·fup ·ζm·ζsr where ·dr is direct solar radiation (W /m2); ·Aiis i-th rowcollectors area (32 m2); ·cos κiis optical losses factor for i-row; ·fup is receptor unlit part factor; ·ζmis mirror nominal reflectivity(0.92); ·ζsr is secondary reflector nominal reflectivity(0.77). The total real solar radiation that affects the receptor is the sum of the real reflected solar radiation of all rows, multiplied bythe total shadowfactor ftot: radrealtot = 11 X j=1 radrealj·ftot. The theoretical solar radiation (W/m2)can beeasily computed without taking into accountthe optical loss factor cos κand the total shadowfactor ftot,and hence radtheotot =dr ·Atot ·ζm·ζsr where Atot is the total reflectivesurface, equal to 352 m2. 4. THERMAL MODEL The thermal model of the plantis focused on heat exchanges in the receptor. All geometrical andstructural data (likelength, exterior diameter, thickness) of absorber pipe, glass cover and secondary reflector areknown as well as thermal parametersnecessary to calculate heat flow, absorptivityand emissivity.In the thermal model, the quantities that havebeen considered as input are: •Ambienttemperature(oC); •Inlet water temperature (oC); •Water flow(m3/s); •Real solar radiation, obtained with the optical model (W/ m2). Fig. 4. Thermal model scheme. The thermal model allows to knowtemperature distribution in the absorption pipeand in heat-transfer fluid (water) along the pipeat agiven moment, as well as the temporary variation of the temperature at determined points of the pipe. For the model the same hypotheses made in Camacho et al. [1997] havebeen adopted. Under these assumptions, byapplying theconservation of the energy in the metal pipeof alength volume dx over atime interval dt,the 18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011 14787 energy balance for eachinfinitesimal volume is obtained. Time variation of internal energy on eachvolume can be expressed as: ∂U ∂t=ρmCmAm ∂Tm ∂t∂x where Cmis the specific heat capacity of the metal,Amis the metal section area,Tmis the temperatureof the metal in Celsius and ρmis the metal density. This variation of internal energy is equals to the difference between the energy of solar radiation incidenton absorber pipeandthe sum of the energy ceded bythe pipeto the ambientand the energy transferred to the fluid. The energy that affects the absorber pipeis the energy due to real solar radiation reflected bycollectors. The expression for solar radiation energy is: Erad =ωr·B·dx where ωris thereal solar radiation and Bis thetotal collector aperture. The energy transferred byabsorber pipeto the ambient and the energy transferred to the fluid haveasimilar expression: Eamb =Hl(Tm−Ta)Bdx Ef=Ht(Tm−Tf)Ldx where Taand Tfare the temperatures of the ambient and of the fluid in Celsius degrees, Lis the inner diameter of absorber pipeand Hl,Htare the coefficient of global thermal losses and the coefficient of metal-fluid transmission, respectively,both expressed in W/m2C.Therefore, the variation of internal energy is equal to: ∂U ∂t=ωr·B·dx −Hl(Tm−Ta)·B·dx −Ht(Tm−Tf)·L·dx The global coefficientof thermal losses is considered per unit of mirror area, that is the power in Watts lost per mirror square meter and Celsius degree. The coefficientof metal-fluid transmission corresponds to the thermal power lost through thecontact surface between fluid and metal. The energy balance in the metal can thus beexpressed as: ρmCmAm ∂Tm ∂t=ωrB−HlB(Tm−Ta) −LHt(Tm−Tf). In asimilar manner, the first lawof thermodynamics can beapplied to the volume of fluid considered: ∂U ∂t=HtL(Tm−Tf)dx −˙m(Hx+dx −Hx) where ˙mis the infinitesimal fluid elementmassvariation, Htis the coefficientof metal-fluid transmission just described while (Hx+dx −Hx)is the differencein enthalpy on the fluid elementof base areaAfand length dx.The enthalpy variation is equal to ∂H ∂xdx =Cf ∂Tf ∂xdx under the assumption that specific heat capacityat constantpressure do not vary in the fluid elementconsidered. Thus, internal energy variation in the infinitesimal element can beindicated as: ∂U ∂t=HtL(Tm−Tf)dx −˙qρfCf ∂Tf ∂xdx where ˙qis the waterflowin m3/s.On the other hand, the internal energy variation in the fluid elementcan be expressed as afunction of the fluid temperature: ∂U ∂t=ρfCfAfdx∂Tf ∂t where Afis the fluid section area, (Af·dx)is the fluid elementvolume, ρfthe fluid densityand Cfthe fluid specific heat capacity.Using this relation,energy balance in the fluid elementhas the following expression: ρfCfAf ∂Tf ∂t+ρfCf˙q∂Tf ∂x=HtL(Tm−Tf). The equations for the energy balance in the infinitesimal elementare therefore the following, where the subscripts mand fare referred respectively to metal and fluid: ρmCmAm ∂Tm ∂t=ωrB−HlB(Tm−Ta)+ −HtL(Tm−Tf) ρfCfAf ∂Tf ∂t+ρfCf˙q∂Tf ∂x=HtL(Tm−Tf). These equations are only applicable to activezones,that is, the parts of absorber pipethat receivethe beam radiation. The equations that describethe performance in apassive elementare similar, except that solar energy inputis nil. Thus, the model for the complete field is built byaseries of activeand passiveelements The equations derived can besolved using an iterative process with finite differences. The temperature of the fluid and of the absorber pipeare calculated for eachtime interval and for eachelement. For the plantconsidered, avalue of 1meter for dx length andatime integration interval of 0.5 seconds havebeen chosen. Atwostage algorithm has been implemented to solve temperature equations. In the first stage the temperatures of the fluid and of the metal are calculated supposing that the fluid is in asteady state. In the second stage the fluid temperature is corrected in function of the net energy transported bythe fluid. 1stStage Tm(n, k)=Tm(n, k−1) +∆t ρmCmAm · ·(ωrB−HlB(Tm(n, k−1) −Ta(k)) −HtL(Tm(n, k−1) −Tf(n, k−1)) Tfp(n, k)=Tfp(n, k−1) + HtL∆t ρfCfAf (Tm(n, k−1) −Tf(n, k−1)) 2ndStage Tf(n, k)=Tfp(n, k) 18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011 14788 −˙q∆t Af∆xTfp(n, k)−Tfp(n−1,k)) In these difference equations, Tf(n, k)and Tm(n, k)are the temperatures of fluid and metal in the segmentn during the k-th time interval, while Ta(k)is theambient temperature during the k-th time interval. Moreover Tfp is the fluid partial temperature before the correction while Tfis the final fluidtemperature. Thus, the outputs of the model are: •Metal temperature of eachpipesegment (oC); •Partial fluid temperature of eachpipesegment (oC); •Final fluid temperature of eachpipesegment (oC). 4.1 Model Identification Although all variables are neededfor acomplete description of the process, the only output of the thermal model that can beused for acomparison to real data is the final temperature of the last segmentof the absorber pipe. In fact there are only twotemperature sensors, one at the beginning and one at the end of absorber pipe, and so is not possible to compare the temperatures of the model internal segments with real ones, neither the metal temperatures of the pipesegments can becompared to the real ones, becausethere are not sensors to monitor them. All geometrical and structural data of the plantare known and thus the constantvalues of the model can bedetermined. The only magnitudes that are unknown are the global losses coefficientHland the metal-fluid transmission coefficientHt.These coefficients are typical of the particular process, and so are valid only for the dynamics of the plantconsidered.Hland Hthavebeen expressed as polynomial functions of second degree of the temperature. The differentparameters of thermal model equations have been determined using real data from the plant. The best fit in the least-squaressense minimizes the sum of squared residuals, aresidual being the difference between an observed value and the value provided byamodel. In this case, the model developed is linear in the parameters and anumerical algorithm for linear model is used to solveit. The model can beexpressed in matrix form as y=Aθ +ǫ, where θis the vector of parameters, yis the vector of measurements, Aasuitable matrix, while ǫis the vector of measurementerrors. 5. VALIDATION AND EXPERIMENTAL RESULTS The model has been validated with real data from the plant. Figures in this section shows the wateroutput temperature given bythe model and the real one for a clear day(Figure 5) and for aparticular windy and cloudy day(Figure 6). As it can beseen, the model prediction curvefollows very well the real one in the clear day,whilein the second case the presence of the wind and of banks of clouds changes the heat transmission condition, issuesthat are not taken into accountin the model. In fact, the wind especially if strong, influences the thermal conductivityto the ambient, increasingthe thermal losses. Moreover the presence of banks of clouds hinders further temperature increase and resultsin aseries of peaks that are still followed well enough bythe model with amaximum difference betweenreal and model data of about8oC. Furthermore, direct solar radiation was estimated from the global radiation measurements and the estimation errors are bigger with clouds. In anycase, the dynamic model developed reproduces the dynamic of the real plantwith enough accuracy to beused for testing control strategies. 11:00 12:00 13:00 14:00 15:00 16:00 17:00 18:00 110 120 130 140 150 160 170 180 Time Water Temperature (°C) Water Outlet Temperature (°C) for the day 27 May 2010 Model Prediction Real Data Fig. 5. Comparison between model predictions and real data of water outlet temperature for May27, 2010. 11.00 12.00 13.00 14.00 15.00 16.00 17.00 18.00 90 100 110 120 130 140 150 160 Time Water Temperature (°C) Water Outlet Temperature (°C) for the day 17 November 2009 Model Prediction Real Data Fig. 6. Comparison between model predictions and real data of water outlet temperature for Nov.17, 2009. REFERENCES A. Larbi, M. Godin, J. Lucas (2000). Analysis of two models of (3D) Fresnel collectors operating in the fixedaperturemode with atracking absorber.Solar Energy, 69, No. 1:1-14. J. Faco and A.C. Oliveira (2009). Numerical simulation of alinear Fresnel solar collector concentrator.SET2009 -8th International Conference on Sustainable Energy Technologies, Aachen (D). E.F. Camacho, M. Berenguel and F.R. Rubio (1997). AdvancedControl of Solar Plants.Springer. M. Spoladore (Supervisors: E.F. Camacho andM.E. Valcher) (2010). Analysis and Modeling of aThermosolar Plant.Automation Engineering Master Thesis. P.Bermejo, F.J. Pino and F. Rosa (2010). Solar absorption cooling plantin Seville.Solar Energy,84: 1503-1512. 18th IFAC World Congress (IFAC'11) Milano (Italy) August 28 - September 2, 2011 14789