A Binary Integer Programming Method for Optimal Wind Turbines Allocation

: The present study introduces a Binary Integer Programming (BIP) method to minimize the number of wind turbines needed to be installed in a wind farm. The locations of wind turbines are selected in a virtual grid which is constructed considering a minimum distance between the wind turbines to avoid the wake effect. Additional equality constraints are also included to the proposed formulation to prohibit or enforce the installation of wind turbines placement at speciﬁc locations of the wind farmland. Moreover, a microscopic wind turbine placement considering the local air density is studied. To verify the efﬁciency of this proposal, a square site was subdivided into 25 square cells providing a virtual grid with 36 candidate placement locations. Moreover, a virtual grid with 121 vertices related with a Greek island is also tested. All simulations conducted considering the area of geographical territory, the length of wind turbine blades, as well as the capacity of each turbine.


Introduction
The modern power grids are complex networks involving a range of energy sources. Traditionally, the production of electric energy took place in the so-called conventional power stations where fossil fuels were used. More recently, the nuclear stations were used in order to increase the amount of produced power energy [1]. However, in the last decade, there was a tendency to adopt renewable energy sources, like wind, sun, and biomass, which are characterized by limited environmental impact. These plants are gradually acquiring significance in many areas and countries. Of all these renewable sources, wind presents the highest penetration and is the most favorable renewable generation option.
Wind turbine generates electricity by making use of renewable wind energy, capturing the wind's kinetic energy and transforming it into electric power. A wind turbine array installed at a distinct geographical area, producing a desired amount of power, is called a wind farm [2]. If the values of the total power required from the wind farm and the capacity of each wind turbine are known, we can determine the maximum number of wind turbines needed to be installed in a specific area. The amount of power generated in a wind farm is dependent on various conditions, such as the morphology of the land, the direction and speed of the wind, and the dimension and characteristics of turbines.
In the last years, the optimal wind turbines allocation problem, known also as the wind farm location optimization problem (WFLO), has been introduced in the literature and several methods have been proposed to solve it [3][4][5][6][7][8][9][10][11][12]. In [3], the optimal locations of wind turbines are calculated using a genetic search code in conjunction with wake superposition. Considering horizontal-axis wind turbines, the 3-D methodology presented in [4] determines their interference effects and performance as well as their flow-field. The interface, enhancement and evaluation of wake and boundary layer models were the main targets examined by the Efficient Development of Offshore Windfarms (ENDOW) project introduced in [5]. The genetic algorithm used in [6] provides the optimal solution of the problem ensuring the minimum number of wind turbines selected and limiting the wind farm territory occupied. In works [7][8][9], the wind turbine placement has been optimally addressed using a simulated annealing algorithm, a genetic algorithm, and a greedy algorithm, respectively. In [10], Monte Carlo simulations are used to optimally arrange the wind turbines in a wind farm. Real-coded genetic algorithms are adopted to solve the optimal placement of wind turbines in [11], while a particle swarm optimization model is introduced for the optimal micro-placement of wind farms in [12]. In [13], the problem is solved considering three cost components of wind turbines dependent on the machine, the tower height, and the generator location. A comprehensive study of WFLO considering both cost and maintenance of wind turbines, the wind speed and direction based on a Weibull distribution, and the wake effect is presented in [14]. An integer nonlinear formulation of the problem is suggested in [15]. Proposal [16] formulates the problem based on a no speeds and coefficients particle swarm optimization (NSC-PSO) method, while in [17] a Firefly Algorithm (FA) is adopted to get the optimal solution. In [18], offshore wind turbines are optimally installed considering the control of the energy efficiency of the wind farm. A probabilistic exhaustive search method considering the Distribution Network Operator (DNO) acquisition market environment allowing wind turbines developers the maximization of the investment's Net Present Value (NPV) over a planning horizon is proposed in [19]. In [20], a meta-heuristic formulation based on the Harmony Search algorithm is proposed to find the optimal locations of wind turbines considering the production efficiency and the capital cost of the investment. Jensen's wake effect model is incorporated to a Mixed Integer Linear Programming (MILP) formulation to find the optimal locations of wind turbines [21]. In [22], the minimum number, locations, and capacities of wind turbines, to be installed in a distribution system, are determined using a commercial package for optimal power flow. In [23], the problem solved using a Particle Swarm Optimization (PSO) methodology along with an optimal power flow algorithm and considering security constraints.
In the proposed method, we subdivide a given predefined territory to squares considering the wake effect and constructing a virtual grid to find the optimum allocated wind turbines for this grid. The avoidance of the wake effect is achieved considering a minimum distance between the wind turbines. The effect of local air density is also considered for the microscopic placement of wind turbines. The problem is formulated as a binary integer programming (BIP) method with binary decision variables. The method can be easily implemented and is successfully tested in a 36-nodes virtual grid as well as in a real case for Skyros island using the MATLAB environment and the easily available bintprog solver [24]. It must be noted that the proposed methodology is quite preferable since the theory is richer, the mathematics are nicer, while the computation is simpler considering linear problems [25].
The main contributions of the proposed wind turbine placement method are: • It considers the wake effect of wind turbines and the local air density.

•
The solution of the problem derived in an efficient and systematic manner.

•
The convergence of the algorithm is ensured using existing commercial optimization software.

•
The execution time, for the cases used for simulations, is reasonable although the problem belongs to the category of planning problems. • It permits the usage of easily available solvers.
The rest of the paper is organized as follows. Section 2 reviews the binary integer linear programming. In Section 3 the proposed method is explained. Section 4 contains numerical results and discussion of performance, followed by concluding remarks in Section 5.

Binary Integer Programming
Linear programming is a mathematical modeling technique used to solve a large variety of problems with moderate effort. A linear program is formulated using a linear objective function of the continuous unknows subject to linear equality and linear inequality Clean Technol. 2021, 3 464 constraints. A binary integer programming problem is formulated as a linear program where the variables are integers. The general binary integer programming formulation can be stated as: where x = (x 1 , x 2 , . . . , x n ) T is the n-dimensional vector of unknowns, x ∈ {0, 1} n , f (x) is the objective function of the problem, h(x) =0 and g(x) ≤0 are vector function equality and inequality constraints, respectively, while S is a subset of n-dimensional space.

Problem Formulation
Considering the energy production from the wind, wind turbines are used. Wind turbines are characterized by the wake effect where cones (wakes) of slower and more turbulent air are created behind them. To avoid this phenomenon, wind turbines must be installed considering a minimum distance between them [14]. If the wind turbines to be installed are constructed by the same vendor and the rotor diameter for each one of them is D ij , the wind farmland can be divided into square cells, having, each one of them, an edge of αD ij meters; α is an integer number greater than one. The binary decision (placement) variables for the problem are defined as The BIP formulation for the optimal allocation of wind turbines can be written as where x is a binary decision variable vector, f (x) is a vector function, whose entries are non-zero if a vertex and its adjacent vertices are candidates for wind turbine placement; 0, otherwise;1 is a vector whose entries are all equal to one; w i is the wind turbine installation cost, and n is the vertices of the virtual grid. In order to show how the elements of vector function f (x) are formulated, we assume a part of the virtual grid as shown in Figure 1. The constraint f i , corresponding to vertex i , is written as follows:

Wind Turbine Placement Prohibition at Specific Locations
Since the grid is virtually constructed, it is possible to have some sites that are not preferable for the installation of wind turbines. The prohibition of wind turbine placement The constraint f i , corresponding to vertex i, is written as follows:

Wind Turbine Placement Prohibition at Specific Locations
Since the grid is virtually constructed, it is possible to have some sites that are not preferable for the installation of wind turbines. The prohibition of wind turbine placement in specific vertices of the grid is related with several parameters such as the different surface conditions, the sheltering effects due to buildings and other obstacles located near the wind farm and the wind variation and intensity which are dependent on by different topographic conditions, e.g., hills and ridges, open plain, sheltered terrain, open sea, and seacoast. Other factors that may prohibit the wind turbine placement is the noise produced by them, the visual impact concerning the neighbouring to private properties or public views, and the influence on areas that have been characterized as archaeological sites or present environmental interest.
When a vertex i is not preferable for the wind turbine installation, the problem of Equation (3) is reformulated as follows: where the equality constraint x i = 0 expresses the prohibition of the wind turbine installation at vertex i, while P is the set of the prohibited vertices.

Mandatory Wind Turbine Placement at Specific Locations
The geographical topology presents an important role for the optimal wind turbine installation. Generally, the sites of placement are characterized by high wind potential. Some vertices of the virtual grid, corresponding to specific locations, may have higher wind potential than others. As a result, they are preferable to wind turbine installation. Another classical case for mandatory placement of wind turbines is the case of complex wind farm terrain, where some sites are adequate for placement, and others are sheltered. The distance from loads to be fed is another factor that affects the wind turbine placement, since the cost of wiring becomes prohibitive as the distance increases. Moreover, in cases of older wind farms, some wind turbines must be relocated inside the wind farm to reduce the financial losses.
The mandatory installation of wind turbines is incorporated into the initial problem formulation (3) as min∑ where the equality constraint x i = 1 expresses the mandatory installation of the wind turbine at vertex i, while L is the set of the preferable vertices.

Microscopic Wind Turbine Placement Considering Local Air Density
For the proper design of a wind farm, many characteristics are considered. One of the most important characteristics is the local climate in terms of the wind speed. Considering the IEC61400-12-1 standard, the normalized wind speed,V n , will be where ρ 0 is the reference air density, typically equal to 1.225 kg/m 3 , and ρ is the local air density at wind turbine hub height. The problem, incorporating the local air density, can be written as min∑ where V n,i and ρ i are the normalized wind speed and the local air density at wind turbine hub height at vertex i, respectively.

Simulation Results
Let an area of 5.0625 km 2 . If the blades of wind turbines are 45 m long and the capacity of each turbine is 1.8 MW, a virtual grid with 36 candidate vertices is constructed [26]. Note that the integer number, α, is selected equal to five (α = 5) [27]. The proposed BIP formulation is tested on a 3.4 GHz Intel(R) Core(TM) i7-2600 processor with 16 GB of RAM using the bintprog solver of the MATLAB optimization tool [24]. The bintprog solver uses a linear programming (LP)-based branch-and-bound algorithm. The solution is derived by solving a series of LP-relaxation problems, in which the binary integer requirement on the variables is replaced by the weaker constraint 0 ≤ x ≤ 1 [24].

Wind Turbine Placement without Considering Equality Constraints
Considering the above-described grid and setting the wind turbine installation cost equal to 1 for every vertex i, w i = 1, ∀i, the proposed BIP formulation (3) is written as follows: The inequality constraints for vertices 6 and 9 are used for illustration: The solution of the above problem gives a total number of ten wind turbines located at vertices 2, 5,9,13,18,22,26,32,35, and 36. The total installed capacity was 18 MW, while the total time needed to get the final result was 0.132 s. Figure 2 depicts the abovementioned solution.
To examine the impact of the integer number α on the performance of the proposed formulation, we conducted simulations considering different values of α for different virtual grid sizes. Table 1 Table 1 includes the ratio of the minimum number of wind turbines to the total number of virtual grid vertices. It is obvious that almost 1/3 of the virtual grid vertices are used for wind turbine installation. Moreover, Figure 3 shows a plot of the objective function value versus iterations to highlight the convergence performance of the proposed formulation for α = 5.
The solution of the above problem gives a total number of ten wind turbines located at vertices 2,5,9,13,18,22,26,32,35, and 36. The total installed capacity was 18MW, while the total time needed to get the final result was 0.132 s. Figure 2 depicts the above-mentioned solution. To examine the impact of the integer number α on the performance of the proposed formulation, we conducted simulations considering different values of α for different virtual grid sizes. Table 1 Table 1 includes the ratio of the minimum number of wind turbines to the total number of virtual grid vertices. It is obvious that almost 1/3 of the virtual grid vertices are used for wind turbine installation. Moreover, Figure 3 shows a plot of the objective function value versus iterations to highlight the convergence performance of the proposed formulation for α = 5.

Wind Turbine Placement Prohibition at Specific Locations
The case of wind turbine prohibition is simulated selecting the vertices 18, 19, and 32 of the virtual grid as the locations in which the wind turbine installation is not permitted. According to (5), we can write:

Wind Turbine Placement Prohibition at Specific Locations
The case of wind turbine prohibition is simulated selecting the vertices 18, 19, and 32 of the virtual grid as the locations in which the wind turbine installation is not permitted. According to (5), we can write: min∑ Clean Technol. 2021, 3

468
The solution of the above problem gives a total number of ten wind turbines located at vertices 3, 6, 7, 11, 16, 20, 24, 27, 31, and 35, while the total installed capacity was 18 MW. The total time needed to get the final result was 0.399 s. Figure 4 presents the problem solution considering the prohibition of wind turbine installation at specific locations. It must be noted that the addition of the three equality constraints in the initial problem formulation increases the time needed to solve the problem and changes the location of wind turbines, while their minimum number remains the same.

Mandatory Wind Turbine Placement at Specific Locations
The exploitation of higher wind density at specific locations in terms of mandatory wind turbine placement is studied using the model described in (6). Considering that vertices 1, 15, and 34 of the grid, as described in Section 3.2, are more preferable for the wind turbine placement, we can write The optimal set includes ten wind turbines located at vertices 1, 4, 7, 12, 15, 23, 25, 26, 30, and 34, presented in Figure 5. The total time needed to get the final result was 0.416 s, while the total installed capacity was 18MW. Once again, the equality constraints added change the wind turbine locations and increase the time needed to solve the problem. It must be noted that the addition of the three equality constraints in the initial problem formulation increases the time needed to solve the problem and changes the location of wind turbines, while their minimum number remains the same.

Mandatory Wind Turbine Placement at Specific Locations
The exploitation of higher wind density at specific locations in terms of mandatory wind turbine placement is studied using the model described in (6). Considering that vertices 1, 15, and 34 of the grid, as described in Section 3.2, are more preferable for the wind turbine placement, we can write The optimal set includes ten wind turbines located at vertices 1,4,7,12,15,23,25,26,30, and 34, presented in Figure 5. The total time needed to get the final result was 0.416 s, while the total installed capacity was 18 MW. Once again, the equality constraints added change the wind turbine locations and increase the time needed to solve the problem. Clean Technol. 2021, 3, FOR PEER REVIEW 9 Figure 5. Optimal solution of the wind turbine placement problem considering mandatory wind turbine placement at specific locations.
The total time needed to get the final result was 0.423 s. To consider the geographic topology and the wind data using a real case, we selected a geographical territory in Skyros island, Greece. Skyros island belongs to Northern Sporades and covers an area of 308km 2 . In this case, we incorporated also the mandatory placement and prohibition of wind turbine's placement. Figure 6 presents in the left hand the wind data as they received from Centre of Renewable Energy Sources & Saving (CRES) [28], while the right hand of the figure presents the geographic morphology of the island. Considering the morphology of the area, the wind potential, the length of wind turbines' blades (each blade is 45mlong), and that the integer number α is equal to 5, a squared area of 20.25km 2 is constructed including 121 candidate vertices (yellow squared area in Figure  6). Considering also that six candidate vertices are excluded due to lower wind potential If the wind speed for each area is 11.00 m/s, 11.20 m/s, 11.10 m/s, and 11.05 m/s, respectively, while the corresponding local air density at wind turbine hub height of each area is 1.199 kg/m 3 , 1.200 kg/m 3 , 1.260 kg/m 3 , and 1.228 kg/m 3 , the normalized wind speed is calculated for each area, respectively. The optimal wind turbine placement problem considering the local air density at wind turbine hub height will be written as The solution of the above problem gives a total number of ten wind turbines located at vertices 3,5,7,16,18,20,25,29,30, and 33, while the total installed capacity was 18 MW. The total time needed to get the final result was 0.423 s.
To consider the geographic topology and the wind data using a real case, we selected a geographical territory in Skyros island, Greece. Skyros island belongs to Northern Sporades and covers an area of 308 km 2 . In this case, we incorporated also the mandatory placement and prohibition of wind turbine's placement. Figure 6 presents in the left hand the wind data as they received from Centre of Renewable Energy Sources & Saving (CRES) [28], while the right hand of the figure presents the geographic morphology of the island. Considering the morphology of the area, the wind potential, the length of wind turbines' blades (each blade is 45mlong), and that the integer number α is equal to 5, a squared area of 20.25 km 2 is constructed including 121 candidate vertices (yellow squared area in Figure 6). Considering also that six candidate vertices are excluded due to lower wind potential and four other vertices are mandatory selected due to the highest wind potential, the implementation of the proposed algorithm gave as a result 42 wind turbines, while the execution time was 1.475 s.
Clean Technol. 2021, 3, FOR PEER REVIEW 10 and four other vertices are mandatory selected due to the highest wind potential, the implementation of the proposed algorithm gave as a result 42 wind turbines, while the execution time was 1.475 s. Figure 6. Wind data and geographic topology of Skyros island, Greece.

Comparison with Other Works in the Literature
In order to demonstrate the efficiency of the proposed algorithm regarding wind turbine placement, Table 2 compares the results for cases described in Sections 4.1 and 4.2, in terms of wind turbines number, locations, and execution time, with the corresponding results of [26]. For the sake of simplicity, we name the case in which equality constraints are not considered in the problem formulation as Case 1 and the case in which the wind turbine placement at specific locations is prohibited as Case 2, respectively. As can be seen, for all cases, the proposed method provides equal number of wind turbines, but in different locations. Moreover, the execution time is quite similar. The difference between the execution times included in Table 1, for Cases 1 and 2, is related with the introduction of additional equality constraints to the problem formulation, increasing the complexity of the problem as well as the corresponding execution times. Given that the rated power of each wind turbine is 1.8 W Μ and considering that the capacity factor is 0.26, the annual energy production of the wind turbine used, will be 0.26·1.8·8760 = 4100MWh/year. Since the optimal solution, for each simulated case, includes 10 identical wind turbines of the same vendor, the annual energy production of the

Comparison with Other Works in the Literature
In order to demonstrate the efficiency of the proposed algorithm regarding wind turbine placement, Table 2 compares the results for cases described in Sections 4.1 and 4.2, in terms of wind turbines number, locations, and execution time, with the corresponding results of [26]. For the sake of simplicity, we name the case in which equality constraints are not considered in the problem formulation as Case 1 and the case in which the wind turbine placement at specific locations is prohibited as Case 2, respectively. As can be seen, for all cases, the proposed method provides equal number of wind turbines, but in different locations. Moreover, the execution time is quite similar. The difference between the execution times included in Table 1, for Cases 1 and 2, is related with the introduction of additional equality constraints to the problem formulation, increasing the complexity of the problem as well as the corresponding execution times. Given that the rated power of each wind turbine is 1.8 MW and considering that the capacity factor is 0.26, the annual energy production of the wind turbine used, will be 0.26·1.8·8760 = 4100 MWh/year. Since the optimal solution, for each simulated case, includes 10 identical wind turbines of the same vendor, the annual energy production of the wind farm will be 41,000 MWh/year. Moreover, considering that the total installation cost, including turbine (ex-works), foundation, electric installation, grid connection, control systems, consultancy, land, financial costs, and road, is 1000 EUR/MW [29], the total installation cost for wind turbines will be 10·1.8·1000 = 18,000 EUR. Table 3 provides a summary of the variables and settings used for the simulations, while a qualitative comparison of the proposed method with other works published in the literature in terms of the contribution is presented in Table 4.  Table 4. Qualitative comparison of the proposed BIP method with other works in the literature.

Method Contribution
The problem is addressed in discrete space by considering a gridded version of the wind farm site and designating the corresponding cell centers as candidate wind turbine locations. The problem is solved using a genetic algorithm.
[4] A scheme for rare positioning of the wind turbines that significantly reduces the land requirements through higher density of turbines is proposed. [5] The primary goal of this work is to link boundary-layer and turbine wake models for the better determination of the wind shear and turbulence profiles inside large offshore wind farms. [6] A genetic algorithm approach is proposed to minimize the number of turbines installed and maximize the production capacity considering non-uniform wind with variable direction, uniform wind with variable direction, and unidirectional uniform wind. [8] A Weibull distribution is used to commonly approximate the wind speed probability density function at a site to correctly calculate the wind turbine power generating capability. [7,9] Different algorithms for the wind turbines placement, namely, gradient search algorithm (GSA), greedy heuristic (GHA), genetic (GA), simulated annealing (SAA), and pattern search (PSA), are examined in terms of the computation time needed to get the optimal solution and the quality of the solution.
[10] An approach based on Monte Carlo simulation to maximize the energy production and minimize installation cost criteria is proposed. [11] This proposal optimizes the locations of each turbine within its cell boundaries in terms of maximizizng the generated wind power considering real-coded genetic algorithms. [12] An approach for wind farm design, based on evolutionary algorithms and related techniques emphasizing particle swarm optimization.
Proposed A given territory is subdivided into squares considering the wake effect and constructing a virtual grid to find the optimum allocated wind turbines for this grid.
The avoidance of the wake effect is achieved considering a minimum distance between the wind turbines. The effect of local air density is also considered for the microscopic placement of wind turbines.

Conclusions
This paper proposes an efficient method for wind turbine placement in power plants, by using a binary integer programming approach. Additional equality constraints are also included to the proposed formulation in order to prohibit or enforce the installation of wind turbines placement at specific locations of the wind farmland. Data related with the area of geographical territory, the length of wind turbine blades, and the capacity of each turbine were used for more-realistic simulations. A microscopic wind turbine placement considering the local air density was also studied. The suggested formulation was tested in MATLAB environment using the bintprog solver and a 36-nodes virtual grid, while a real case scenario concerning the Greek island of Skyros was also examined. The simulation results were compared in terms of wind turbines number, locations, and execution time, with those of a previous author's work. The results showed that the proposed methodology can be effectively implemented considering or prohibiting the installation of wind turbines in specific locations. The execution time needed is comparable with that of a previously published work, while the minimum number of wind turbines is the same.