A Numerical Methodology for Simple Optimization of Marine Hydrokinetic Turbine Array

Tidal stream energy, due to its high level of consistency and predictability, is one of the feasible and promising type of renewable energy for future development and investment. Applicability of Blade Element Momentum (BEM) method for modeling the interaction of turbines in tidal arrays has been proven in many studies. Apart from its well-known capabilities, yet there is scarcity of research using BEM for the modeling of tidal stream energy farms considering full scale rotors. In this paper, a real geographical site for developing a tidal farm in the southern coasts of Iran is selected. Then, a numerical methodology is validated and calibrated for the selected farm by analyzing array of turbines. A linear equation is proposed to calculate tidal power of marine hydrokinetic turbines. This methodology narrows down the wide range of turbine array configurations, reduces the cost of optimization and focuses on estimating best turbine arrangements in a limited number of positions.


Introduction
Along with increasing popularity of renewable energy resources in the 21st century, as a viable alternative for fossil fuels, research and development in the field of marine renewable energies gained a significant additional boost. Among different types of marine renewable energies, tidal energy, due to its high level of consistency and predictability, is one of the most feasible and promising sources of energy for future development and investment. Tidal energy resources may be exploited potentially (i.e., tidal potential energy), or kinematically (i.e., tidal stream energy). Focus of this study is on tidal stream energy devices and their performance in a farm.
One important step towards commercial deployment of tidal energy is performing feasibility studies. There are several approaches to mapping and characterization of tidal stream energy resources in a real geographical site. Uihlein and Magagna reviewed 14 mapping studies across different sites of the world [1]. One numerical approach is to use 2D model built based on MIKE 21 to simulate tidal hydrodynamics. The approach was used by Gao et al. [2], Zheng et al. [3] and Zhang et al. [4] for tidal current resource assessment in China coats. Also, Radfar et. al. [5] utilized the technique for mapping of tidal stream energy resources in the southern coasts of Iran. Their are used as input of this study.
Significant harvestable kinetic energy via an array of Marine Hydrokinetic (MHK) turbines in a tidal channel accelerated its utilization in global energy sector. Hence, performing extensive investigation on fluid mechanics of tidal energy farms gained significant attention during recent years. Array layout should be optimized such that there is a balance between maximising total power capture and power per turbine while minimising impact on natural flow (i.e., do not slow the natural flow by a significant fraction) [6]. Studies on finding an optimized layout of tidal turbine array directly influences economic benefits of tidal farms. Complexity and range of possible array layouts leave a vast optimization space that could be explored. Therefore, proposing new numerical approaches to decrease the need for extra simulations is of great importance for future studies in this field and especially if they could effectively decrease computational cost. It is the main goal of this study.
Applicability of Blade Element Momentum (BEM) method, as an appropriate numerical methodology for this type of investigations was proved in many studies. Historically, BEM was developed by Zori and Rajagopalan [7] for its application in modeling flow fields around helicopter rotors and investigating its performance. BEM simulates aerodynamic effect of rotating blades via time-averaged body forces over a full revolution of blades, placed on a thin full continues disk with an area equal to the swept area of the blades [8]. Many researches suggest to couple BEM with Reynolds Averaged Navier-Stokes (RANS) equations, in which RANS equations are employed to model the flow field within the entire computational domain. Masters et al. [9] compared numerical modeling techniques for tidal stream turbine analysis including Blade Element Momentum Theory (BEMT), Actuator Disk model, RANS-BEM (or BEM-CFD), blade-resolved moving reference frame and coastal models based on the shallow water governing equations. They showed that BEM-CFD is a promising tool to investigate the near-field, wake characteristics, recovery and interaction. In this approach, BEM method is used to model the turbine, and CFD is employed to model flow properties elsewhere in the domain.
There are many parameters that have impact on turbines in tidal farms, such as lateral distance between turbines (lateral layout), inter-row spacing (longitudinal layout), rotation direction of rotors, and so on. Malki et al. [10] investigated the influence of tidal stream turbine spacing on performance of group of three tidal stream turbines using BEM-CFD model and found that an optimum downstream spacing of the third turbine depend on the across-stream lateral spacing between two upstream turbines. Turnock et al. [11] did research on wake modeling and analysis of tidal turbine arrays to improve the efficiency of energy harvesting process using a coupled RANS-BEMT approach. They concluded that for a fixed number of turbines minimising the lateral spacing within each row, with a small number of staggered rows spaced as longitudinally as far apart as practical, is the most effective strategy for energy capture. There is a limit for minimisation of this lateral spacing and a minimum tip-to-tip distance must be between two adjacent turbines in a row. Malki et al. [10] using BEM-CFD stated that for a lateral spacing of 3 diameters or greater, there was little difference in downstream recovery, however, for smaller spacing, recovery was much slower. They also observed that in staggered layout, for lateral spacing of 2 diameters or greater, flow accelerates between upstream turbines and the performance of the downstream turbines is maximised by reducing the downstream distance to benefit from such flow acceleration.
Harison et al. [12] investigated The effect of ambient turbulence on the inter-row spacing and reported that an increased ambient upstream turbulence decreases the wake recovery length and thus the required inter-row spacing decreases. Also, if it is practical, it is proposed that the downstream distance increases to 75% of wake velocity recovery distance to ensure consistent wake recovery and thus consistent power extraction of downstream devices [13,14]. From practical point view, it is not beneficial to increase downstream distance unlimitedly, because reducing the distance between turbines has several advantages shorten the submarine cable length, efficient utilization of ocean space, lower the difficulty of tidal power farm construction, and so on [15]. This clearly signifies the importance of optimisation studies using computationally efficient numerical models. Mozafri [16] proposed general numerical methodology based on fitting a linear relationship between available Kinetic energy flux 2R upstream the turbine and the extracted power in the absence of the real geometry of the MHK turbines. The present study aims to generalize this step-by-step methodology and to resolve in conceptual problem in its proposed linear relationship.
This research explores a fundamental linear relationship to estimate the power capture of full scale turbines using available kinetic energy flux. In this regard, a real site for developing a tidal farm in the southern coasts of Iran is selected. Then, the numerical methodology of Mozafri [16] is being applied and modified for the selected farm. Arrays of turbines are analyzed and a linear equation developed for estimating the tidal power of devices. The novelty of current research is to consider full scale turbines with the presence of real geometry of rotors and presenting a simple linear method to estimate tidal energy output. The article is structured as follows: in Section 2 discussions on tidal site selection are being presented. Detailed description of the numerical model, BEM model and the characteristics of the computational domain are presented in Section 3. In Section 4, the numerical methodology is being calibrated for the investigation of selected tidal stream farm in the selected site. An step-by-step straightforward numerical methodology based on fitting a linear relationship will be presented in this section. Finally, the summary of the conclusions of the research is presented in Section 5.

Tidal Site Selection
In this work a tidal site is selected in the southern waters of Iran. The geometry and the flow field of this site is idealized and numerically modeled to examine its potential for implementation of the first pilot tidal stream energy farm in Iran. As previously noted, for efficient utilization of ocean space, shortening the submarine cable length and lowering the difficulty of tidal power farm construction it is preferable to build a tidal farm with minimum employed area. Besides, There are also other important parameters for the selection of a permitted tidal farm in the prescribed region (i.e. southern waters of Iran): 1. Availability of giant oil and gas resources in this part of Iran and the need for significant capital investment for the tidal energy projects; small pilot farms are preferred. 2. Existence of fragile marine life and habitat in this site; leads to permitting issues for development of large tidal farms in this site. Again small pilot farms are preferred. 3. Minor variations of bathymetry in this site [5]. This will significantly decrease computational time and cost. Furthermore, it smoothen the pathway for idealization of the geometry and the flow field within the site. According to these discussions and based on previous study of Radfar et al.  Table 1.  As can be seen, the selected tidal site geometry is idealized to have the shape of a rectangular cube with the previously mentioned dimensions (see section 2) and a constant depth as shown in Figure 2. It is more straightforward to present farm dimensions based on the diameter of rotors.

Numerical Model
In the present study, the fluid mechanics of the channel is modeled by solving RANS equation closed by k- turbulent model. Operation of turbines within the modeled farm is simulated via BEM coupled with RANS. In this section of paper, a brief description of numerical modeling is provided.

Turbine Modeling: BEM
As mentioned before, BEM is used to numerically simulate operation of turbines in the tidal farm. With BEM turbine blade span is divided into small elements from root to tip. Here, elemental lift and drag forces are computed based on local flow field velocity, computed as solution to RANS equations, and provided lookup table of lift and drag coefficients as a function of Angle of Attack (AOA) for each element. Dynamic effect of blades rotation is simulated based on time-averaged body forces over a thin full disk with an area equal to the swept area of the turbine, over a full revolution of blades [8]. Final time-averaged body forces has the same magnitude, but different sign to the lift and drag forces acting on each segment. This exclusion of the actual geometry of blades from numerical modeling of the turbine blades reduces computational costs even up to 100 times, while providing an accurate effect of the turbine far wake and turbine-turbine wake interaction [9,[18][19][20][21]. This makes the BEM a promising numerical method for the current work.

Boundary Conditions
The imposed boundary conditions are uniform streamwise velocity with magnitude of 1.75 m/s (this velocity is based on the real tidal stream velocity of the selected farm based on the previous study [5]) and atmospheric pressure at inlet and outlet respectively. The other boundaries of computational domain are modeled as slipfree walls. BEM discretization has about 3.86  10 6 elements. The number of nodes along the length, width and depth of the domain are 338, 127 and 29 respectively. Mesh is structured in most of the computational domain except for the regions close to the nacelles. For the sake of achieving mesh independency, the mesh resolution near the rotor zone and nacelles changed in seven steps to get desired mesh accuracy in power estimation of rotors. Attentions should be paid for the meshing of rotor zone. Rotor zone is a thin volume just upstream of the rotor face. Duplicate faces must be avoided, here. Besides, the mesh must have a full 360 unit circular face as the rotor face.

Numerical Results
Before discussing on the results of numerical simulations, equation for Available Kinetic Energy flux (AKE) should be defined. AKE can be expressed as [5]: where  is density of seawater (here, 1025 kg/m 3 ) and V is the tidal current speed [m/s]. AKE is calculated based on the area-averaged velocity in a circular surface with the same area that of rotor face: In this section, discussion on the results is based on Normalized Available Kinetic Energy flux (NAKE), which is calculated as: where V  is the uniform streamwise velocity entering the computational domain through inlet (here, 1.75 m/s). The two other parameters are global and local efficiencies. The global efficiency (E G ) is defined as the power extracted by the turbine in a certain array position, normalized by the power available in the undisturbed flow at the inlet of the channel: The important aspect of the efficiency in the array optimization process is how much kinetic energy flux is going through a given position (array placement efficiency) and how efficiently each device in the array extracts energy from the available kinetic energy at that position (the turbine "local efficiency"). The local efficiency (E L ) of is calculated as: where V local the velocity averaged on the projection of the rotor swept area on a plane one diameter upstream of the turbine.
Finding optimum location for placement of tidal turbines in a farm causes more power output and less depreciation cost for the devices. Although there were several attempts for layout optimization of lab scaled tidal farms, but still steps are to be passed to reach a level where it is possible to implement optimization studies in full scale. However, it is known that strong hardwares with several cores and huge computational cost is needed in this regard. Having these in the mind, current study is very early step toward the goal of having optimized full scale tidal farm. As it is mentioned previously, Mozafari [16] developed a general methodology based on the results of sixty different simulations using various turbine configurations. This methodology supports the idea that there is a linear relationship between the power extracted by the turbine and AKE flux at upstream the device. It was found that all data falls along the same slope which is equivalent to the local efficiency of all turbines being constant, regardless of their position in the array. Owing to this fact, a new methodology has been proposed in the following section to put aside the dependence on the site specifications and generalized it for various engineering applications in different areas all over the world where potential tidal stream power is considerable. In the following, the methodology for the turbine array optimization process is presented step by step: First step: Due to natural (farm width) and technical limitations (in this study, the tip-to-tip distance is 1.5R), maximum of three turbines can be placed at each row. Figure 3 shows the normalized streamwise velocity contour in the middle plane (z=0). It can be concluded that the streamwise velocity increases nearly two times between two adjacent turbines, due to positive interactions of rotating blades. Also, flow accelerates around the middle turbine and thus its captured power is 2% more than the other two. The extracted power, the kinetic energy at 1D upstream, the local and global efficiencies of all three turbines are summarized in Table  2.  Second step: As mentioned earlier, due to the presence of the linear relationship between available kinetic energy at 1D upstream and the extracted power, a methodology had been developed for array optimization objectives. Previous methodology suffers one major problem; they fitted "Y=mX+b" to data, where Y is the extracted power and X is the available kinetic energy flux at 1D upstream. Suppose that X=0; then you can see using this equation Y=b, which is physically meaningless. Generalizing this methodology to account for all natural conditions and removing site requirements from the provisions, modifications have been made and the new methodology is proposed in the current study: "In BEM, it is possible to fit a linear relationship in the form of Y=mX between the extracted power (Y) and the available kinetic energy flux at 1D upstream (X), where m is the mean local efficiency of first row turbines." Owing to this methodology, a linear relationship to estimate the power of downstream turbines based on the results of Table 2 can be calculated, as where, A is the rotor disk swept area and V is the local averaged velocity. Flowchart of numerical methodology is presented in Figure 4. Third step: Now, it is necessary to verify the validity of proposed methodology. In this regard, a second row must be added in order to compare simulated power of turbines against the estimated values of proposed linear relationship. In this step, a suitable location for the turbines of the second row must be specified. In choosing this location paying attention to three issues is of great importance: 1. Basically, it is preferred that the turbines located in a position with maximum kinetic energy flux and thus maximum energy. 2. Turbulence effect should not to be too destructive to cause damages on rotor blades and other devices. 3. A minimum distance between rotors should be adopted as to let the biggest marine biota pass safely through rotors. 4. Downstream turbines cannot be placed very close to the upstream turbines due to the negative effects of rotating blades. These effects can deteriorate tidal devices and eclipse the efficiency of the turbine. At first glance, it seems that the best location for the second row is right after the nozzles of the first row. But, this location could not be reported definitely without further investigations on turbulence effect. In this regard, several questions remain unanswered at present. Further research should be undertaken to investigate the destructive effects of the turbulence forces on the routine operation of elements of the turbines. In general, there should be a minimum distance to avoid harmful effects of wake on downstream devices. To study the performance of the downstream wakes, NAKE flux and Normalised Turbulence Kinetic Energy (NTKE) in different downstream distances are reported in Table 3. It is worth mentioning that NAKE corresponds to the amount of the extracted power and NTKE describes the intensity of the turbulence forces.  Table 3 show that moving from 2D to 5D downstream, available kinetic energy flux and turbulence kinetic energy decrease 7.2 and 16.5 percent, respectively. These findings suggest that there are uncertainties in the selection of the optimized location for the downstream rows. A trade-off should be made between low turbulence effects and high expected power.
Hitherto, the first and second conditions from the three mentioned conditions for the selection of the location of the downstream row have been described. Regarding to the third issue, the biggest marine mammal in the farm location must be detected. For the selected farm location in this study, the biggest marine mammal is Black dolphin, for safety passage of this biota, a minimum distance of six meters should be considered between two adjacent turbines. This helps reducing threatening environmental risks. In Figure 5, seven different routes for placing the corners of downstream turbines are determined. Looking at the contours of streamwise velocity clearly shows that staggered positions (L1, L3, L5 and L7) are desirable due to bigger tidal current velocity and hence, more AKE. Three of these four candidate positions should be selected for the placement of the second row turbines. According to the first simulation, the rotor 13 has more power than rotor 11, turbines of the second row are oriented to the left side of the farm. Thus, the corners of the second row turbines should be placed on the passes L1, L3 and L5. Besides, to avoid disruptive effects of rotating blades, a minimum distance of 1.5R has been applied between ends of the nozzles of the first row with the beginning of the nozzles of the second row. Figure 6 shows the velocity contours of this simulation. Also, the quantitative results are summarized in Table 4. Now a comparison must be made between the estimated power based on Eq. (1) and the actual power based on the simulation. Results of the seventh column of Table  4 clearly shows that that the maximum difference does not exceed 3%.  More validation tests must be performed for a solid conclusion on the viability of the proposed formula of the methodology. To validate the proposed linear relationship (Eq. (2)) further simulations have been conducted in 21 different locations downstream of the first row turbines. All these results are summarized in Table 5. As it is seen here, the maximum difference is about 6%, which is acceptable and desirable in engineering applications. An interesting observation here is that the extracted power from staggered turbines is significantly more than coaxial cases. This clearly signifies our previous conclusions. The same approach can be applied for the placement of any other downstream rows.

Discussions and Conclusions
This study set out to develop a numerical methodology based on the linear relationship to reduce the computational cost and to obviate the need for mesh and modeling procedures of the downstream rows. BEM is used to numerically simulate operation of turbines due to its computational efficiency. In tidal turbines like any other structures and devices, cyclic load can cause fatigue damage and structural cracks in long-term. Turbulence is defined as an additional fluctuating velocity component added to the mean velocity. This fluctuating velocity acts as a dynamic load on the turbine with many cycles within a second of time [22]. DNV guideline states that changes in the local flow field and in the ambient turbulence intensity are to be expected due to the wakes of upstream turbines. These wake effects imply a considerably increased turbulence. This effect may be significant even when the spacing between the tidal turbines in the tidal array is as large as 10 rotor diameters (longitudinally). It is reported that the wake effects can even dominate fatigue effects [23].
With current knowledge and lack of commercial data on structural effects, there will be uncertainties in selecting the best position for downstream turbines. The findings of this research suggest that in this regard, several questions remain unanswered at present. Further research should be undertaken to investigate the destructive effects of the turbulence forces on the routine operation of elements of the turbines. Since, yet there is no commercial scale tidal farm in the world, for better understanding an example of existing wind farm is presented here. Due to same nature of HAHTs and wind turbines, this example helps better understanding of the trade-off. Figure 7 shows the Lillgrund wind farm in Sweden. As can be seen, 3.3D lateral distance and 4.3D longitudinal distance have been implemented to avoid destructive effects of wakes. Therefore, it is expected that such a distance should be applied in tidal farms to avoid compact layout. [24].

Figure 7: Layout of Lillgrund wind farm in Sweden
In Section 3, a linear relationship between AKE flux at 1D upstream and the expected power of the downstream turbines was fitted and a step-by-step straightforward numerical methodology has been proposed. This methodology narrows down the wide range of turbine array configurations, reduces the cost of optimization and focuses on estimating best turbine arrangements in a limited number of rows and columns. Decreasing the need for extra simulations is of great importance for any optimization studies in this field and significantly reduces the computational cost.
The main novelty of this linear methodology is its undeniable advantage in engineering applications. A simple linear relationship, reduces the efforts for meshing and the numerical modeling of the downstream turbines. Whereof, there will be always a growing need for more number of turbines to extract the potential tidal energy in the farms, the remarkable value of the methodology in industrial estimations will be better apparent. Especially, the results of this research showed that the difference between the predicted value and the actual power does not exceed 6%. So, this numerical methodology is capable of attracting great attention in practical applications.