Solar Central Receiver Systems, SCRS in what follows, are an interesting and promising technology in the field of solar-based renewable energies. They were initially proposed in 1957 and, since then, they have been under continuous development (Gordon, 2013; p. 618). This kind of power facilities is especially characterized by its scalability, operating efficiency and output stability (Avila-Marin et al., 2013; Besarati and Goswami, 2014; Collado and Guallar, 2012, 2013).

In simple terms, SCRS consist of a set of high reflectance mirrors, called ‘heliostats’, and a radiation receiver which is usually on top of a tower. Heliostats track the apparent movement of the Sun throughout days to concentrate solar incident radiation over the receiver. Consequently, it is under high-density solar radiation. Then, a heat-interchange process is carried out in it to heat a working fluid. This fluid can be finally used in a classic turbine cycle to generate electric energy. Due to the high temperatures reached, thermodynamic performance is high (Besarati and Goswami, 2014; Collado and Guallar, 2012). In Fig. 1, a simple SCRS facility is depicted. Further information about them can be found in Alexopoulos and Hoffschmidt (2017), Behar et al. (2013), Camacho et al. (2012), Müller-Steinhagen (2013), Reddy et al. (2013), Stine and Geyer (2001).

This work deals with optimizing the heliostat field of SCRS facilities, which is also known as ‘collector subsystem’. The set of heliostats and their configuration define this part of facilities. This subsystem can cause up to 40% of operating energetic loss and can take up to 50% of the initial investment according to Besarati and Goswami (2014), Jones et al. (2007). Therefore, it must be carefully optimized when designing a SCRS power facility. Several optimization objectives are possible. For instance, optical efficiency maximization, as studied in Besarati and Goswami (2014), Noone et al. (2012), is one of them. Production cost reduction, as mentioned in Ramos and Ramos (2012), Sánchez and Romero (2006), is another valid option. Finally, some variations and combinations of the previous ones, as studied in Zhang (2007), are also common objective functions.

As summarized in Carrizosa et al. (2014), the design of heliostat fields is usually based on two basic strategies. The first one is to apply and optimize a certain heuristic distribution pattern. The classic radial staggered model, proposed by the University of Houston (Stine and Geyer, 2001; Chap. 10), and the biomimetic spiral, recently proposed by Noone et al. (2012), are two examples of patterns. In Mutuberria et al. (2015), where several patterns are compared, more examples and references can be found. The second strategy is to select final heliostat positions from a defined set (as iteratively done in Sánchez and Romero, 2006). However, a mixed approach is usually applied: An over-sized field is initially generated by following a particular distribution model to get a set of available positions over the ground. Then, the best ones are finally selected until the real power requirements are fulfilled. In fact, this strategy is quite common and it has been used, for instance, in Besarati and Goswami (2014), Noone et al. (2012), Ramos and Ramos (2012), Pitz-Paal et al. (2011), Sánchez and Romero (2006). This latter approach permits more degrees of freedom to get optimal fields.

In this work, parametric optimization of the biomimetic spiral pattern is studied. The design objective is to maximize the optical efficiency of fields. This distribution scheme could be ultimately applied to both directly generate final fields and in a hybrid design methodology. Thus, it is important to rely on an efficient and reliable optimizer to adjust its design variables. In this context, two stochastic population-based algorithms, called micraGA and UEGO, are applied to solve the problem at hand. Their structure seems to be especially suitable to solve it. However, as far as the authors know, none of them had been previously used in this situation. Hence, the main goal of this work is to study their behaviour in it. Additionally, a Pure-Random Search algorithm is used to study the complexity of the problem and to get a baseline perspective of the quality of solutions. Finally, a Brute-Force Grid is applied to approximate where the optimal solution is deterministically. In contrast to the active trend of deploying large fields with several thousands of heliostats (e.g. 4120 units in Khi Solar One and 173500 ones in Ivanpah), this paper considers small and medium size fields. Specifically, facilities with several hundreds of heliostats, like CESA-I and PS10, with 300 and 624 units respectively, define the scope of interest. This is due to their lower requirements of space and initial investment, i.e. higher applicability. However, the size of the different problems is gradually increased, which makes it possible to discover general patterns. The results obtained seem scalable.

The rest of the paper is organized as follows: In Section 2, the optimization problem and the objective function are exposed. In Section 3, the different optimizers applied to the problem are described. In Section 4, the experiments carried out and the results obtained are shown. Finally, in Section 5, conclusions are drawn and future work is exposed.

The design of heliostat fields, as introduced, will be aimed at maximizing yearly irradiance weighted efficiency, η y e a r , I . It is a common objective function which is used, for instance, in Besarati and Goswami (2014), Noone et al. (2012). It is focused on the optical efficiency of the heliostat field, while also considering the available solar radiation at every instant through the year. This magnitude can be formulated as Noone et al. (2012): (1) η year , I = ∑ day = 1 365 ∫ sunrise sunset I b ( t ) η ( t ) d t ∑ day = 1 365 ∫ sunrise sunset I b ( t ) d t where I b ( t ) is the instantaneous beam irradiance and η ( t ) is the instantaneous optical efficiency. It must be noted that Eq. (1) is defined in abstract terms. Hence, it must be configured by selecting the appropriate model for every encapsulated concept. These aspects are detailed in Noone et al. (2012), where the unweighted version of the equation is also included.

The optical efficiency of a field, η ( t ), measures its performance in redirecting and concentrating incident solar radiation over the receiver at every instant t (Stine and Geyer, 2001; Chap. 10). In this context, an instant is equivalent to a certain solar position and incident radiation density on a particular day and time. It defines the heliostat field configuration. As an efficiency factor, η ( t ) is defined between 0 and 1 (minimum and maximum possible efficiency, respectively). It is internally formed by a selected set of sub-factors which model different sources of energy loss in heliostat fields. Discarding any of them is equivalent to assume it as 1, i.e. not causing energy loss. For the scope of this work, η is defined for every heliostat as in Noone et al. (2012): (2) η = η cos · η s b · η itc · η a a · η ref where η cos , η s b , η itc , η a a and η ref are the cosine, shadowing and blocking, interception, atmospheric attenuation and reflectivity efficiency, respectively. These concepts will be summarized later in this section while also including the particular model applied to their calculation. Thus, η ( t ) in Eq. (1) is referred to the average obtained from the values of every active heliostat. Further information about these efficiency factors can be found in Besarati and Goswami (2014), Stine and Geyer (2001), Noone et al. (2012). Additionally, the alternative definition of Eq. (2) used in Collado and Guallar (2012) is very descriptive. Its nomenclature not only includes time and the coordinates of heliostats but also the effect of their neighbours.

Instantaneous beam solar radiation, I b ( t ), weights pure optical efficiency. As introduced, it adds an energetic perspective to Eq. (1). It modifies the instantaneous optical efficiency of a certain field considering the availability of solar radiation at every instant studied. This factor is also abstractly defined and the selected model will be described later in this section.

The instantaneous component, t, encompasses a year, which is also divided into days. Days are finally considered as raw sets of instants, i.e. apparent solar positions. The way in which days are discretized is also abstract and configurable. For instance, in Noone et al. (2012), some interesting notes about considering a variable time step are presented. However, those strategies are out of the scope of this work and a constant time step will be used to compare the selected optimizers.

Figure 2 gives an overview of the problem to solve. Eq. (1), which directly depends on the heliostat field design, is the objective function to be maximized. It can be seen as simulation procedure of this part of SCRS facilities. A common context defined by the number of time steps per year, the configuration of every internal model, etc., is assumed. Therefore, it is necessary to define the design of candidate fields to be analysed with the underlying objective function while optimizing. As introduced, fields will be designed by directly applying the biomimetic spiral pattern proposed in Noone et al. (2012). Consequently, the optimizers applied will need to find the parameters of this pattern, ‘a’ and ‘b’, that generates a field with a maximum value of Eq. (1). As usual in engineering problems, the objective function does not directly depend on the design variables. Next, the internal models selected for Eq. (1) and the biomimetic pattern are described.

The optimization problem described has been solved with four different algorithms: a genetic algorithm (called micraGA), a memetic global optimizer (called UEGO), a Pure-Random Search (that will be referred as PRS) and a Brute-Force Grid (named as BFG in what follows).

The problem context has been implemented in C++. The objective function, which requires simulating the candidate fields throughout the year, has been parallelized due to its significant computational cost. The simulation of instants has been distributed among concurrent threads as proposed in Cruz et al. (2017). The POSIX-Threads (PThreads) library has been used for parallelization. The optimizers described in Section 3 have also been implemented in C++. The compiler used is g++ 4.8.1 with optimization level ‘O2’. The experiments have been carried out in a cluster node that features an Intel Xeon E5 2650 with 16 cores and 64 GB of shared RAM. The evaluation of the objective function can hence deploy up to 16 concurrent threads.

Regarding the problem definition, the latitude of the facility is 37.4 ∘ North. The height of the tower is 100 m with a cylindrical receiver of 10.5 m of height and 8.5 m in diameter. Since a north-field configuration is required, only those positions of the pattern at the north of the tower are considered. All heliostats are assumed to have a plain reflective surface of 100 m 2 (10 m × 10 m) whose center is at 5 m over the ground. The time step throughout days is constant. Specifically, for any studied field, Eq. (1) is computed at every hour from sunrise to sunset. For the selected solar model and latitude, approximately 4000 instants can be uniformly generated in a year. However, the solar positioning model shows symmetry, which is inherited by radiation density estimation. Thus, simulations only evaluate about 2400 real instants.

The considered problem sizes range from 50 to 500 heliostats to place. According to Besarati and Goswami (2014), a and b, are defined in [ 2.0 , 8.0 ] and [ 0.45 , 0.70 ], respectively. For BFG, the discretization steps in a and b are 0.05 and 0.005, respectively. This configuration has also been used to plot the shape of the objective function shown in Fig. 8 for both the smallest and the biggest instances. As can be seen, the search space is quite homogeneous, which is mainly due to keeping it limited to a robust distribution model. UEGO and micraGA have been configured, after preliminary experimentation, with a robust set of parameters for all the studied problems. They have been configured to get the best results possible, with stability in spite of their stochastic nature, and as fast as possible. For UEGO, this is achieved with a limit of field evaluations equal to 1000, 20 levels of search, a limit of 15 species and a minimum radius equal to 0.0001. micraGA has been configured to have a population size of 70 and to select 70 progenitors from tournaments of 12 individuals at every cycle. micraGA will execute a total of 4 cycles. Its cooling procedure will divide mutation increments by 9 . 0 ( x − 0.5 ) , where x is the cycle (when it is larger than 1). Mutation amplitude in a is in range [ − 1.0 , 1.0 ] while it is in range [ − 0.15 , 0.15 ] in b. Finally, PRS has been adjusted to take a similar time to UEGO and micraGA. Specifically, it is allowed to execute 1250 cycles. This time constraint would be interesting in case of using this process to generate preliminary fields in a more complex design method.

Table 1 contains the results obtained after experimentation. Each row represents a particular problem size while the sets of columns show the results achieved by every optimizer. The records of BFG contain the efficiency achieved and the runtime for every problem size. However, due to their non-deterministic nature, the records of UEGO, micraGA and PRS consist of the average efficiency achieved, the standard deviation and the runtime for every problem size. The experiments with non-deterministic methods have been repeated for 10 times. Additionally, the last row includes the average of all the instances. The best efficiency in every case is in bold font.

As can be seen in Table 1, the achieved efficiencies are numerically very similar for the same instance. Besides, multiple combinations of a and b result in fields of similar efficiency. In Pitz-Paal et al. (2011), where heliostat field optimization is also studied, a similar situation was observed. Generally speaking, the real impact of such variations in efficiency would ultimately depend on the field size. In fact, the differences between those methods with a better average (UEGO and BFG) and the others tend to be higher with more heliostats. If we studied the final average efficiency of every method as that of a virtual field, the difference between the best and worst methods would be of a few thousands of euros a year (less than 6,000 euros). However, if these results were scaled to larger fields, the same degree of variation would represent an extra benefit of several thousands of euros (less than 60,000 euros). These values would not be determinant if compared to the total budget of the facility. Nevertheless, they could be achieved by simply selecting the most appropriate optimizer at design. Moreover, since analytical models might tend to overestimate efficiencies because of their simplifications, it is important to work with high precision.

UEGO gets the best average efficiency in all the studied cases. It is also very stable as can be observed in its values of standard deviation. In fact, it usually returns the same solution for each problem instance even after different executions. They all are also in the region where the optimal solution is expected to be according to the preliminary analysis with BFG. Therefore, UEGO is neither affected by its stochasticity nor by the multimodal structure of the objective function.

BFG achieves efficiencies quite near to those obtained by UEGO. However, its runtime is significantly higher. Hence, its integration in complex design procedures would be less practical. On the one hand, BFG is up to 10 times slower than UEGO with a relatively coarse discretization of the search space while even obtaining slightly worse results. On the other hand, BFG is deterministic and its resolution can be arbitrarily adjusted through its discretization steps. These properties make it theoretically interesting for studies such as the one presented in this work. However, as can be deduced from Table 1, BFG would need an enormous amount of time to handle larger problem instances, wider search spaces and higher accuracy. More sophisticated methods such as UEGO outperforms BFG with less computational effort.

Regarding micraGA, its performance is quite low when compared with UEGO. In fact, it is even marginally outperformed by PRS in two cases, the instances of 200 and 450 heliostats. This behaviour is consistent with the idea that the complexity of micraGA is between those of UEGO and PRS. Furthermore, it is not particularly stable if compared with the standard deviation obtained by PRS and especially by UEGO. micraGA does not converge adequately due to the structure of the search space. It would need a significantly higher runtime to get closer results to UEGO. Thus, since micraGA cannot outperform UEGO and it is complex if compared to PRS, it is not a good option to be further considered.

Finally, despite being the worst ones, the results obtained by PRS are quite interesting. In spite of its simplicity and lack of orientation during the search, depending on the precision needed, it could achieve acceptable results. This is mainly due to two reasons: i) its reduced logic overhead makes it possible to study numerous combinations in a similar runtime to UEGO and micraGA; ii) it is not affected by the structure of the search space. In fact, PRS could be interesting as part of a more complex heliostat field optimization method because of its simplicity. However, its instability and lack of robustness compared to UEGO make it not appropriate as a standalone method.

In this work, the complexity of heliostat fields design has been studied. An optimization problem aimed at maximizing yearly irradiance weighted efficiency has been defined. The design procedure followed directly tries to optimize the parameters of the promising biomimetic distribution pattern. Four optimization algorithms have been studied to solve the problem, namely, BFG, UEGO, PRS and micraGA. The possibility of coupling their use within more complex methods is also discussed considering computational cost. Whereas Genetic Algorithms and BFG are state-of-art methods for the present problem, UEGO had not been applied to solve it before. PRS, as far as the authors know, could also be included in this set. Since the objective function is wrapped within a distribution robust pattern, it is heavily multimodal.

UEGO is the best option out of the considered ones. This method gets the best solutions in a stable and fast way for all the studied instances. Hence, it should be considered for heliostat field design based on the described distribution pattern. It also seems to be fast enough to be included as part of a more complex design procedure. The field efficiencies obtained by BFG are very similar to those achieved by UEGO. In fact, for a given resolution, BFG could obtain the best solution possible due to the extensive search that it performs. However, it would require significantly more runtime than UEGO. Thus, BFG does not seem to be appropriate for large fields, wide search spaces or to be included in a complex procedure of several stages. Regarding PRS, it turns out to be an interesting option considering its simplicity. Although it should not be applied to generate final fields, PRS could be appropriate in a multi-stage procedure depending on the precision requirements. Finally, the genetic optimizer proposed, micraGA, does not seem to be a good option. It does not explore the search space properly. Its results are of low quality and unstable if compared to those of UEGO. In fact, they are very similar to those of PRS in spite of being more sophisticated than that method.

There are some possible expansion points for future work. First, a better alternative to micraGA could be explored. Its design should aim at a trade-off between the complexity of UEGO with the best results and the simplicity of PRS with acceptable quality. Second, it could be interesting to compare the performance of UEGO and PRS in a multi-stage optimization process. Finally, extending this comparison context to encompass several distribution patterns and timing variations in the objective function would be of great interest.