An optimization tool to design the field of a solar power tower plant allowing heliostats of different sizes

An optimization tool to design the field of a Solar Power Tower plant allowing heliostats of different sizes E. Carrizosaa,b, C. Dom´ınguez-Bravoc,∗, E. Fern´andez-Carad,b, M. Queroe aDepartment of Statistics and Operations Research, University of Seville, Spain bIMUS - Mathematics Institute of the University of Seville, Spain cBCAM - Basque Center for Applied Mathematics, Bilbao, Spain dDepartment of Differential Equations and Numerical Analysis, University of Seville, Spain eAbengoa Solar New Technologies, Seville, Spain Abstract The design of a Solar Power Tower plant involves the optimization of the heliostat field layout. Fields are usually designed to have all heliostats of identical size. Although the use of a single heliostat size has been questioned in the literature, there are no tools to design fields with heliostats of several sizes at the same time. In this paper, the problem of optimizing the heliostat field layout of a system with heliostats of different sizes is addressed. We present an optimization tool to design solar plants allowing two heliostat sizes. The methodology is illustrated with a particular example considering different heliostat costs. Keywords: solar thermal power tower, field layout, multi-size-heliostat field, heuristic algorithm, greedy algorithm 1. Introduction Solar power tower (SPT) system is known as one of the most promising technologies for producing solar electricity due to the high temperatures reached that result in high thermodynamic performances; some reviews on solar thermal electricity technology are [1,2,3,4]. In an SPT system, direct solar radiation is reflected and concentrated by the heliostat field in a receiver placed at the top of the tower. At the receiver, the solar energy is converted into thermal energy by heating a fluid which can be used to generate electricity through a conventional thermodynamic cycle. The heliostat field is composed of a group of mirrors having usually two-axis tracking system to reflect the direct light from the sun to the receiver aperture. ∗BCAM - Basque Center for Applied Mathematics, Mazarredo, 14 E48009 Bilbao, Basque Country, Spain. E-mail: [email protected] Preprint submitted to International Journal of Energy Research December 30, 2015 The optimization of the field layout to minimize the levelized cost of energy (LCOE) is a challenging problem due to several reasons: the problem is of very large dimension (with hundreds, or even thousands, of variables), involves nonconvex constraints and the objective function is nonsmooth, hard to compute and multimodal [5]. With the purpose of reducing the complexity of the problem, a geometrical pattern is frequently imposed, i.e. the heliostats locations follow a fixed distribution. Usually, a radially-staggered [6,7,8,9,10], spiral [11] or grid [12] distribution is used. Thus, the heliostat field layout is calculated through the optimization of a low number of parameters defining the geometry of the selected distribution. Since 1970s, different research programs have been financed to study the heliostat design aiming to break down the heliostat costs while maintaining the collected energy. Original ideas have been developed in recent studies (hexagonal [13], bubble [14], minimirror array [15] and other geometries [16,17,18]), see also [19,20,21,22,23]. Some promising prototypes are the following: the autonomous heliostat (CIEMAT PCHA project [24,25]), the EASY heliostat (hEliostat for eAsy and Smart deploYment [26]) and the SCS5 eSolar next generation heliostat [27]. The heliostat optical design influences the overall performance of the system (ratio of ground coverage, number of heliostats, receiver size and tower height), and it is influenced by the cost (manufacturing and assembly processes, canting, installation, calibration, etc.), and wind loads among others, see [14,28,13,29, 18]. As pointed out in several studies [14,21,26], efficient fields can be designed with large or big heliostats (148 m2ATS [14], 121 m2Sanlucar [30], 116 m2 Sener [21]) but also using small or micro heliostats (16 m2AORA Solar and Heliko-DLR [14], 7.5m2[31], 4.3m2SHP-CSIRO [29]). As pointed out in [32], the use identical sizes may not lead to optimal fields. Despite this, the design of fields using heliostats of different sizes together remains, as far as the authors are aware of, unexplored. In this paper we will focus on the optimization of a multi-size-heliostat field (a heliostat field with different heliostat sizes) using a pattern-free method. For simplicity, we will assume that the tower and receiver are fixed (to address the whole optimization problem, see [33]), the pedestal height remains the same for all heliostat sizes and, as usual, all heliostats will be assumed to focus onto the same target point: the aperture centre. The rest of the paper is organized as follows. In Section 2, we describe the main ingredients affecting the behaviour and performance of the SPT system. Our methodology to solve the optimization problem is explained in Section 3. In Section 4, we present the heliostat sizes used in this paper and we apply the proposed algorithm and analysis tools to a typical plant design. Finally, in Section 5, our main results are summarized and some perspectives for further work are presented. 2 2. Problem statement In this section, we explain the meaning of the variables involved in the optimization process. We also present the constraints that have to be satisfied, the cost and energy functions (which are the elements to be considered for the computation of the objective function LCOE), and the optimization problem itself. 2.1. Variables The heliostats locations, given by the coordinates (x, y) of their centres, and the heliostat sizes d, are the variables to be used. From now on, we will denote by Ω the collections of coordinates of the centres and sizes of the heliostats, namely (x, y, d). The set Ω is described as follows: Ω = {(xi, yi, di) for i∈[1, N] with (xi, yi)∈ S and di∈ D} ,(1) where Ndenotes the total number of heliostats, Sis the set of heliostats coordinates and Dis the set of heliostat sizes. We assume that the set Dis finite. For simplicity, all heliostats are assumed to be rectangular and have the same pedestal height, although they can have different dimensions. Note that these assumptions help to reduce the computation of the shading and blocking effects caused by large-size heliostats on the smaller ones. 2.2. Constraints Usually, when designing an SPT system, a fixed time is used to evaluate the plant operation. This time is known in the literature as the design point, denoted here by Td. Let ΠTd(Ω) be the power input obtained at the design point. Then, a minimal power input has to be achieved, that is, the following constraints has to be satisfied: ΠTd(Ω) ≥Π0.(2) Due to technical reasons, the heliostats must be located within a given region S0⊂R2: S ⊂ S0.(3) The heliostats located in the field have to rotate freely avoiding collisions with other heliostats. Consequently, we have to include constraints forcing the heliostats not to overlap: ||(x, y)−(x0, y0)|| ≥ δ(d) + δ(d0)∀(x, y, d),(x0, y0, d0)∈Ω with (x, y)6= (x0, y0),(4) where the radius of the clear-out circle for heliostat size dis δ(d)=0.5diag(d)+ 0.5ds. Here, diag(d) denotes the heliostat diagonal and dsis a positive constant, related to installation errors and heliostat accessibility, which remains equal for all the heliostat sizes in this paper. 3 2.3. Functions The cost and the annual energy are the functions involved in the objective function. The cost function C=C(Ω) takes into account the investment in power plant equipment (tower, receiver and heliostat field), purchasing of land and civil engineering costs: C(Ω) := K+ Ψ(Ω),with Ψ(Ω) := X d∈D c(d)Nd,(5) where Kis a constant including all fixed costs (independent of the configuration of the heliostat field) and Ψ(Ω) represents the heliostat field cost function. The number of heliostats of each size is denoted by Ndand c(d) denotes the cost per heliostat of size d. All costs associated with the heliostats (mirror modules, support structure, drives, pedestal, foundation, field wiring, etc.) are included in c(d) and, for simplicity, they are supposed to be independent of the heliostat position. The annual energy input function E=E(Ω) takes the form: E(Ω) := ZT 0 I(t) N X i=1 ϕ(t, xi, yi, di,Ω) dt , (6) where I(t) is the so-called instantaneous direct solar radiation and ϕrepresents the product of the heliostat efficiencies (usual in this framework), that is, ϕ= fref fat fcos fsb fsp . Specifically, fref is the heliostat reflectance factor, fat is the atmospheric efficiency, [34,35]; fcos is the cosine efficiency, [35]; fsb is the shading and blocking efficiency [6,36], and, finally, fsp is the interception efficiency or spillage factor [37]. The annual energy of the plant is computed with a procedure similar to NSPOC (Nevada Solar Power Optimization Code) [38]. We refer the reader to [34,37,35] for further details. We have developed a Matlab prototype to adapt the energy calculation when having different heliostat sizes and, in particular, address the shading and blocking effects. 2.4. Optimization Problem The optimization problem we are addressing can be written as follows: (P)          Minimize F(Ω) = C(Ω)/E(Ω) Subject to ΠTd(Ω) ≥Π0 Ω⊂ S0× D ||(x, y)−(x0, y0)|| ≥ δ(d) + δ(d0) for (x, y, d),(x0, y0, d0)∈Ω with (x, y)6= (x0, y0). (7) In this problem, the number of heliostats is not fixed in advance. Note that, even fixing this number, the huge amount of heliostats in recent commercial plants makes this problem very difficult to solve, as pointed out in [2]. 4 Some of the heliostat efficiency functions depend on the heliostat area and/or the position in the field (interception efficiency [39,40,5], atmospheric efficiency [2,11], etc.) Hence, the heliostats annual energy per unit area values are different depending on their positions, see Figure 1. These values are similar in general for both sizes although they have a different behaviour in regions below the two quadrant diagonals (due to the interception efficiency, see Figure 2). 3. Field optimization algorithm The aim of this paper is the field layout design when having different heliostat sizes. The proposed algorithm aims to work with any selected size having arbitrary aspect ratio, cost, etc. We describe our approach in the case of two sizes (big and micro heliostats), but the methodology extends easily to the general case. The proposed procedure, called Expansion-Contraction Algorithm, starts with a large-size heliostat field (calculated following the Greedy Algorithm explained in Section 3.1, see also [33,41]) and complements it by inserting smallsize heliostats. Following this algorithm, large-size heliostats will be located at the best positions taking advantage of the most favourable region near the tower. Then, two consecutive phases, called Expansion and Contraction, are applied and repeated until a stopping condition is fulfilled. At the Expansion Phase, small-size heliostats are inserted with the Greedy Algorithm. At the Contraction Phase, the best heliostats are selected according to their LCOE per unit area values and the worst are sequentially deleted. The Expansion-Contraction algorithm is explained in detail in Section 3.2. In order to allow the possibility of mixed-fields when having big and micro heliostats, the algorithm starts locating big heliostats first. If smaller heliostats were considered as first candidates, only infeasible positions would remain for the location of large heliostats if required. This way, the same idea presented in [42] is followed: “ The best strategy to fill a case with stone, pebble and sand is as follows. First filling the case with the stones and then filling the gap left from the stones with pebbles and in the same way, filling the gap left from pebbles with sand. Since filling in opposite direction may leave the stones or pebbles outside.” 3.1. Greedy Algorithm The procedure presented in this paper makes use of the Greedy Algorithm, designed in [33] to calculate a field layout with a single heliostat size. This is a method that sequentially locates the heliostats one by one in the field at the best feasible position. The annual energy values are modified at each step due to the shading and blocking effects that the new heliostat produces in the field. For simplicity, it is assumed that the heliostat cost is independent of the heliostat location, so that the annual energy function Ecan be viewed as the objective function. The heliostats are located freely, without any pre-arranged 5 (West) y coordinate (East) (South) x coordinate (North) Annual Energy Values per unit area 0 0.01 0.02 0.03 0.04 0.05 0.06 Figure 1: Annual energy per heliostat unit area of small-size (GW Hth/m2). Outer semicircle has a radius of 9.95 tower heights (1 km). Annual Interception Efficiency Values 0.2 0.4 0.6 0.8 Figure 2: Interception efficiency (values in [0,1]), overlapping of small-size (thin lines) and large-size (thick lines) heliostat effects. 6 distribution. Only two geometrical constraints have to be taken into account: the field shape constraint (3) and the constraints to avoid heliostat collisions (4). Obviously, the first problem involves locating the first heliostat centre when only the field shape constraint is considered. This problem has an easy-to-handle objective function, because of the absence of shading and blocking effects. In return, when we have already located k−1 heliostats and we have obtained a field Ωk−1={(x1, y1, d1),...,(xk−1, yk−1, dk−1)}that fulfils (3) and (4), the problem Pkdescribed below is difficult to solve, since non-convex constraints are involved and the energy function has a complex shape due to the shading and blocking effects. Let us introduce the notation Ωk= Ωk−1∪ {(x, y, d)}, where (x, y) denotes the variables with respect to which we maximize in problem Pk. Now, we focus on the problem of finding the optimal location of a new heliostat: Pk   Maximize E(Ωk−1∪ {(x, y, d)}) Subject to (x, y)∈ S0 ||(x, y)−(x0, y0)|| ≥ δ(d) + δ(d0)∀(x0, y0, d0)∈Ωk−1. (8) When k > 0, this problem becomes multi-modal due to the shading and blocking effects. A multi-start procedure is used to avoid local minima starting from several randomly selected feasible positions. The final solution is chosen according to the annual energy given by each configuration. The final number of heliostats is given by the algorithm. It stops when the power requirement are reached. 3.2. Expansion-Contraction Algorithm The Expansion-Contraction algorithm starts with a feasible large-size heliostat field that reaches the power input constraint (2) and then makes a series of Expansion-Contraction steps. The Expansion-Contraction algorithm is described in Algorithm 1. The Expansion Phase consists of oversizing the large-size field using smallsize heliostats until a prescribed power input value Π+ 0, greater than Π0, is reached. The small-size heliostats are located one by one following the Greedy Algorithm, recalculating the shading and blocking effects at each step. Small-size heliostats are expected to fill-in possible holes between the large-size heliostats already located due to their smaller area. Moreover, as can be seen in the contour lines shown in Figure 2, they reach higher energy per unit area values in lateral regions. Once the oversized multi-size-heliostat field is obtained, the heliostats are arranged according to their LCOE per unit area values. At the Contraction Phase, the heliostats reaching lowest values are (sequentially) deleted and the number of selected heliostats is determined by (2). This phase has to follow a sequential procedure because once a heliostat is deleted, the shading and blocking effects over its neighbours change and, therefore, their values have to be recalculated. This process can be carried out selecting carefully the active neighbours in order to avoid the recalculation of the annual energy of the whole field and, 7 consequently, reducing the computational time. Oversizing and selection are well-known in the field layout problem, as they are usually used in combination with some fixed-pattern strategies, see [6,7,10,11]. Algorithm 1 Expansion-Contraction Algorithm Require: Π0and Π+ 0 Ω0←Create initial field using large-size heliostats with Greedy Algorithm. Stop when Π0is reached. F0←F(Ω) Υobj ←F0 Υfield ←Ω0 stop ←0 k←0 while k≤kmax &stop = 0 do Expansion Phase: Ω+ k←Oversize Ωkusing small-size heliostats with Greedy Algorithm. Stop when Π+ 0is reached. Contraction Phase: Ω− k←Sort Ω+ kaccording to: LCOE per unit area. Select the best heliostats until Π0is reached. Update: k←k+ 1 Fk←F(Ω− k) Ωk←Ω− k if Fk≥Υobj then Υobj ←Fk Υfield ←Ωk else stop ←1 end if end while return Υfield 4. Results In our experiments we consider the large-size heliostat of the kind HLarge, whose area is 121.34 m2(similar to the usual heliostat used with the selected tower-receiver configuration, see Table 2, the heliostat Sanlucar120 [30,43]). The small-size heliostats are of the HSmall kind, with area 4.35 m2(similar to the SCS5 eSolar heliostat [27]). In Table 1and Figure 3both are fully described. Several studies, see [14,20], support a reduction on the heliostat cost per unit area for small heliostats compared to large heliostats. In this paper, the 8 heliostat cost per unit area is set to 158.61 $/m2and two different cost scenarios are studied: •Scenario-100: For small-size and large-size heliostats, the costs per unit area are identical. •Scenario-80: For small-size heliostats, the cost is only 80%. The corresponding LCOE functions are respectively denoted F100 and F80. For simplicity, in this paper the pedestal height and safe distance remain the same for all heliostat sizes. These assumptions help to reduce the shading and blocking effects caused by large-size heliostats over the smaller ones. However, note that the selected sizes have different aspect ratio, which implies a different clear-out ratio. Therefore, these effects and the computed solutions will depend strongly on the selected heliostat sizes. Heliostat Parameter Large-size Small-size Name HLarge HSmall Width [m] 12.84 3.21 Height [m] 9.45 1.36 Optical height z0[m] 5.17 5.17 Diagonal [m] 15.94 3.49 Safe distance ds[m] 1.70 1.70 Security distance δ(d) [m] 17.64 5.19 Aspect Ratio (width/height) 1.36 2.36 Total Area Ad[m2] 121.34 4.35 Relative Area 1 0.25 Hel. Cost Scenario-100 [$/m2] 158.61 158.61 Hel. Cost Scenario-80 [$/m2] 158.61 126.89 Table 1: Parameter Values (Heliostat sizes) In view of the tower parameters detailed in Table 2, we see that the maximum thermal energy value is reached at coordinates (74.41,0), i.e., 0.74 tower heights to the North. Note that this value belongs to the interval [0.5,1], given in [8]. Therefore, as can also be appreciated in Figure 1, this region (near the tower) is in principle the most favourable to locate heliostats and a higher density of heliostats is expected there, as pointed out in [11,29]. The annual energy per unit area generated by one single heliostat is very similar for both sizes in general. However, the interception efficiency values differ, specially in regions below the two quadrant diagonals, and furnish better results with the small-size heliostat (thick lines in Figure 2). The Expansion-Contraction Algorithm described in Section 3.2 has been implemented in Matlab c , using the fmincon routine to solve the involved optimization subproblems. The specific values for the receiver-field parameters are shown in Table 2. The power input required at the design point Π0is set to 45.03 MWth. The value for the upper limit Π+ 0is set to 49.51 MWth (an increase of 10% 9 [12] M. S´anchez and M. Romero. Methodology for generation of heliostat field layout in central receiver systems based on yearly normalized energy surfaces. Solar Energy, 80(7):861–874, 2006. [13] P. Schramek and D. R. Mills. Heliostats for maximum ground coverage. Energy, 29:701–713, 2004. [14] J. G. Kolb, S. A. Jones, M. W. Donnelly, D. Gorman, R. Thomas, R. Davenport, and R. Lumia. Heliostat Cost reduction study. Technical Report SAND2007-3293, Sandia National Labs., 2007. URL http: //prod.sandia.gov/techlib/access-control.cgi/2007/073293.pdf. [15] J. G¨ottsche, B. Hoffschmidt, S. Schmitz, M. Sauerborn, R. Buck, E. Teufel, K. Badst¨ubner, D. Ifland, and C. Rebholz. Solar Concentrating Systems using Small Mirror Arrays. Solar Energy Engineering, 132:011003–011007, 2010. [16] W. Landman. Sensitivity analysis of a curved heliostat profile. In Proceedings of Annual Student Symposium 2012 in CRSES, 2012. [17] T. R. Mancini. Catalog of solar heliostats. Technical Report III– 1/00, SolarPaces, 2000. URL www.fika.org/jb/resources/Heliostat% 20Catalog.pdf. [18] C. Zang, Z. Wang, H. Liu, and Y. Ruan. Experimental wind load model for heliostats. Applied Energy, 93:444–448, 2012. [19] K. R. Bhargav, F. Gross, and P. Schramek. Life cycle cost optimized heliostat size for power towers. Energy Procedia, 49:40–49, 2014. SolarPaces 2013. [20] J. B. Blackmon. Parametric determination of heliostat minimum cost per unit area. Solar Energy, 97:342–349, 2013. [21] J. Coventry and J. Pye. Heliostat cost reduction – where to now? Energy Procedia, 49:60–70, 2014. SolarPaces 2013. [22] L. Meng, Z. You, A. F. M. Arif, and S. Dubowsky. Shape optimized heliostats using a tailored stiffness approach. Journal of Solar Energy Engineering, 136(021017):(9 Pages), 2014. [23] A. Pfahl. Survey of heliostat concepts for cost reduction. Journal of Solar Energy Engineering, 136(014501):1–9, 2014. [24] G. Garc´ıa Navajas, A. Egea Gea, and M. Romero. Performance evaluation of the first solar tower operating with autonomous heliostats: PCHA project. In Proceedings of the SolarPaces 2004, 2004. [25] A. Pfahl, M. Randt, C. Holze, and S. Untersch´ ’utz. Autonomous lightweight heliostat with rim drives. Solar Energy, 92:230–240, 2013. 16 [26] A. Monreal, M. Burisch, M. Sanchez, D. P´erez, C. Villasante, E. Olabarrieta, D. Olasolo, and A. Olarra. EASY: An innovative design for cost effective heliostats/solar fields. Energy Procedia, 49:174–183, 2014. SolarPaces 2013. [27] P. Ricklin, M. Slack, D. Rogers, and R. Huibregtse. Commercial readiness of esolar next generation heliostat. Energy Procedia, 49:201–208, 2014. SolarPaces 2013. [28] A. Pfahl, M. Buselmeier, and M. Zaschke. Wind loads on heliostats and photovoltaic trackers of various aspect ratios. Solar Energy, 85:2185–2201, 2011. [29] P. Schramek, D. R. Mills, W. Stein, and P. Le Livre. Design of the Heliostat Field of the CSIRO Solar Tower. Journal of Solar Energy Engineering, 131 (024505):(6 Pages), 2009. [30] Solucar. 10 mw solar thermal power plant for southern spain. Technical report, Solucar, Inabensa, Fichtner, CIEMAT and DLR, 2006. URL http: //www.trec-uk.org.uk/reports/ps10_final_report.pdf. [31] M. Renzi, C. M. Bartolini, M. Santolini, and A. Arteconi. Efficiency assessment for a small heliostat solar concentration plant. International Journal of Energy Research, 39:265–278, 2015. [32] L. Crespo, F. Ramos, and F. Mart´ınez. Questions and answers on solar central receiver plant design by NSPOC. In Proceedings of SolarPaces 2011, 2011. [33] E. Carrizosa, C. Dom´ınguez-Bravo, E. Fern´andez-Cara, and M. Quero. A heuristic method for simultaneous tower and pattern-free field optimization on solar power systems. Computers & Operations Research, 57:109–122, 2015. [34] F. Biggs and C.N. Vittitoe. The HELIOS model for the optical behavior of reflecting solar concentrators. Technical Report SAND76–0347, Sandia National Labs., 1976. URL http://prod.sandia.gov/techlib/ access-control.cgi/1976/760347.pdf. [35] F. J. Collado and J.A. Tur´egano. Calculation of the annual thermal energy supplied by a defined heliostat field. Solar Energy, 42:149–165, 1989. [36] G. Sassi. Some notes on shadow and blockage effects. Solar Energy, 31(3): 331–333, 1983. [37] F. J. Collado and J.A. Tur´egano. An analytic function for the flux density due to sunlight reflected from a heliostat. Solar Energy, 37:215–234, 1986. [38] L. Crespo and F. Ramos. NSPOC: A New Powerful Tool for Heliostat Field Layout and Receiver Geometry Optimizations. In Proceedings of SolarPaces 2009, 2009. 17 [39] F. J. Collado. One point fitting of the flux density produced by a heliostat. Solar Energy, 84:673–684, 2010. [40] F. W. Lipps. Four different views of the Heliostat Flux Density Integral. Solar Energy, 18:555–560, 1976. [41] E. Carrizosa, C. Dom´ınguez-Bravo, E. Fern´andez-Cara, and M. Quero. Optimization of multiple receivers solar power tower systems. Energy, 90: 2085–2093, 2015. [42] M. Ergen, C. Coleri, and P. Varaiya. Qos aware adaptative resource allocation techniques for fair scheduling in OFDMA based broadband wireless access systems. IEEE Transactions on Broadcasting, 49:362–370, 2003. [43] Fichtner. Assessment of technology options for development of concentrating solar power in south africa for the world bank. In The Climate Investment Funds, 2010. URL https://www-cif. climateinvestmentfunds.org/sites/default/files/Presentation% 20-%20WB%20(Eskom)%20Project%20-%202010_12_07%20.pdf. [44] R. Osuna, V. Fern´andez, S. Romero, M. Romero, and M. S´anchez. PS10: a 11.0-MWe Solar Tower Power Plant with Satured Steam Receiver. In Proceedings of SolarPaces 2004, 2004. [45] R. Buck. Heliostat field layout improvement by nonrestricted refinement. Journal of Solar Energy Engineering, 136(021014):(6 Pages), 2014. [46] P. Schramek and D. R. Mills. Multi-tower solar array. Solar Energy, 75(3): 249–260, 2003. [47] L. Alon, G. Ravikovich, M. Mandelbrod, U. Eilat, Z. Schop, and D. Tamari. Computer-based management of mirror-washing in utility-scale solar thermal plants. In Proceedings of ASME 2014 (8th International Conference on Energy Sustainability and 12th International Conference on Fuel Cell Science, Engineering and Technology), volume 1, July 2014. 18