Performance of non-uniform tidal turbine arrays in uniform flow

Theoretical models suggest that in order to maximise their collective power output, tidal turbines should be arranged in a single cross-stream row and optimally spaced to exploit local blockage effects. However, because it is assumed that the turbines within these arrays are identical, such models do not consider the possibility of enhanced power production through the exploitation of spanwise variations in local blockage and resistance. In this paper, we use depth-averaged numerical simulations to investigate whether the performance of a tidal turbine array can be further enhanced by varying solely the local blockage, solely the local resistance, or both local blockage and resistance together, across the array width. Our results suggest that for an initially uniform flow field, the optimal tidal turbine array is also uniform, that is to say that it comprises turbines of equal size, spacing, and resistance. This finding is encouraging because it is more cost-effective and much simpler to design each turbine to be the same and to operate in the same way. Together with earlier findings, these results also suggest a more general, and perhaps unsurprising, conclusion that tidal turbine arrays perform best when designed to match site-specific natural flow conditions.


Introduction
In designing a wind or tidal turbine, one of the key challenges is to understand how the performance of the turbine is defined by its interactions with the flow field. This problem is complicated by the need to describe fluid-structure interactions over multiple length scales; a task which becomes exceedingly complex for turbines within large arrays Vennell et al. 2015). One way to approach this problem is to use a simple theoretical model to analyse the performance of an idealised turbine in an idealised flow field. The simplest of these models is actuator disc theory, which approximates the spinning turbine rotor as a porous disc of uniform resistance (Burton et al. 2001). Simplifying the turbine structure in this way eliminates the need to resolve flow features at smaller scales and allows arguments of mass, momentum, and energy to describe how factors such as the geometric blockage ratio and resistance of the turbine affect its performance.
One of the earliest applications of actuator disc theory produced what is commonly termed the Betz limit; a well-known theoretical result which establishes that a wind turbine can extract no more than ∼ 59.3% of the kinetic energy of the air which would pass through the swept area in the absence of the turbine (Betz 1920;Joukowsky 1920;Okulov and van Kuik 2012). Despite the simplifications involved in its development, the Betz limit has proven a useful benchmark for wind turbine performance and has inspired numerous extensions to the classical turbine model. In recent years, actuator disc theory has provided valuable insights into the performance of tidal stream turbines, having been extended to approximate more closely the flow conditions that such turbines experience (Garrett and Cummins 2007;Whelan et al. 2009;Houlsby et al. 2008;Vennell 2010;Draper et al. 2016) and to analyse their performance when placed in different arrangements Willden 2012, 2013;Vogel et al. 2016;Draper and Nishino 2014a, b). Vennell (2010) has shown that in order to maximise their collective power output, tidal turbines should be placed sideby-side to create a single cross-stream row. By maximising the global blockage ratio, which is defined as the ratio of array swept area to channel cross-sectional area, this arrangement allows the turbines not only to extract more power (Garrett and Cummins 2007) but to do so more efficiently, with less power lost in wake mixing (Draper et al. 2010). By making a number of additional assumptions, Willden (2012, 2013) and Vogel et al. (2016) have shown that the performance of a cross-stream row can be further enhanced by adjusting the lateral spacing between the turbines to optimise the local blockage ratio, which is defined as the ratio of turbine swept area to local flow cross-sectional area. Further analysis by Draper and Nishino (2014a, b) suggests that a single, carefully spaced, cross-stream row is in fact the optimal arrangement for tidal turbines, ensuring both a higher power output per turbine and approximately equal distributions of thrust and power among the turbines. The two-scale actuator disc model of Nishino and Willden (2012) appears, therefore, to provide the simplest description of an optimal tidal turbine arrangement. However, because the two-scale model assumes that each turbine within the array is identical and that the depth of the idealised channel is uniform across its width, it does not consider the possibility that more power could be produced by exploiting the spanwise variations in local blockage and resistance which may naturally arise in practice.
In this paper, numerical simulations are used to investigate whether the performance of a tidal turbine array can be further enhanced by exploiting such variations. The combined array and channel model developed by Bonar (2017) is first used to simulate steady, uniform, and depth-averaged flow through an idealised channel with low background roughness. A subgrid-scale actuator disc model is used to introduce an array of global blockage 0.1 and the turbines are arranged to achieve near-optimal uniform local blockage and resistance profiles. Non-uniform profiles are then used to investigate whether the performance of the array can be further enhanced by varying solely the local blockage, solely the local resistance, or both local blockage and resistance together, across the array width. This analysis extends the works of Hunter et al. (2015) and Adcock (2015) to consider the effects of spanwise variations in local blockage as well as local resistance, and to measure array performance in terms of both the collective power output of the turbines and their global power coefficient. The analysis also complements the work of Draper et al. (2016) by examining the performance of a turbine array idealised by a non-uniform local resistance in a uniform flow field, as compared to that of a turbine represented by a uniform local resistance in a non-uniform flow field. Bonar (2017) has developed an idealised numerical model to explore the potential for local blockage effects to enhance the performance of turbines in tidal channels. In the present study, the same model is used to analyse the performance of non-uniform tidal turbine arrays in uniform flow. The main features of the model are described as follows.

Channel-scale flow
Channel-scale flow is simulated by using the open-source hydrodynamic model ADCIRC to solve the depth-averaged shallow water equations by means of a discontinuous Galerkin (DG) finite element scheme (Kubatko et al. 2006(Kubatko et al. , 2009. As illustrated in Fig. 1, the domain of interest is an idealised channel of depth 20 m, width 4 km, and length 20 km. The channel walls are set to allow tangential slip and a steady current is produced by establishing a fixed head difference ζ between the two ocean boundaries, which are placed in deep water (in this case, in water 1 km deep) and positioned far upstream and downstream of the channel (at distances of 50 km) in order to minimise reflections (Adcock 2015). The contours in Fig. 1a describe a typical variation in free surface elevation, as flow is driven from the upstream ocean boundary at x/d = − 3000, through the central channel located between x/d = − 500 and x/d = + 500 (where the depth drops slightly below the still water level), and toward the downstream ocean boundary at x/d = + 3000. The contours in Fig. 1b show the corresponding variation in depth-averaged velocity magnitude, with continuity ensuring that velocities are highest in the shallow central channel and practically negligible at the ocean boundaries where the depth is much greater. Figure 1b also demonstrates the formation of an array-scale wake around the single row of turbines located at x/d = 0, and shows that this wake extends beyond the end of the shallow channel and into deeper water.
In the absence of turbines, the flow loses energy to seabed drag, changes in cross-section, and turbulent mixing. The drag due to seabed roughness is calculated as F = ρ A b u|u|C d , in which ρ is the fluid density, A b is the plan area of the seabed, u is the depth-averaged velocity vector, and C d is a dimensionless seabed drag coefficient. For large-scale tidal models such as these, C d values of ∼ 0.0025 are typical (Soulsby 1997) but, in this study, a much lower seabed drag coefficient is chosen to ensure that the results from the numerical model are comparable to those of the frictionless two-scale actuator disc theory. A value of C d = 0.0005 is found by iteration to provide a good compromise between minimising the effects of channel background roughness and maintaining model stability. A relatively small head difference of ζ = 0.05 m, which produces an unexploited channel velocity of u 0 ≈ 0.964 m/s, is then selected to ensure that the model remains stable even for very high turbine resistance. Mixing is controlled by a spatially and temporally constant horizontal eddy viscosity coefficient, the value of which is calculated, following Borthwick and Barber (1992) and Kuipers and Vreugdenhil (1973), Fig. 1 Plan views of the model domain showing (near) steady-state contours of: a surface elevation η(x, y) normalised by channel still water depth d and b depth-averaged velocity magnitude u(x, y) normalised by unexploited channel velocity u 0 . A steady current is driven from left to right by a fixed head difference ζ and the turbine array is extended inward from one side of the channel (colour online) total depth of flow, with d the still water depth and η the free surface elevation above still water level.

Local-scale flow
The extraction of energy and resulting changes to the flow field are simulated using the open channel actuator disc model derived by Houlsby et al. (2008). Following Draper et al. (2010) and Draper (2011), Serhadlıoglu (2014) introduced this actuator disc model into the DG-ADCIRC code at sub-grid scale, thereby enabling idealised turbines, each defined by a local blockage ratio and local resistance coefficient, to be inserted between numerical elements within the computational grid. The channel-scale flow between these elements provides boundary conditions for the turbine model within their shared edge, which calculates the extracted power and imposes the associated loss of momentum as a discontinuous reduction in fluid depth (Draper et al. 2010). Coupling the DG-ADCIRC and actuator disc models in this way allows the largely two-dimensional array-scale flow problem to be solved numerically, whilst the highly three-dimensional local-scale problem is modelled analytically at sub-grid scale. Though simplistic, this line sink modelling approach has been shown by laboratory experiment to be a reasonably accurate means of describing the momentum deficit imparted by model-scale rows of porous discs (Draper et al. 2013). The local blockage ratio is defined, following Nishino and Willden (2012) and Vogel et al. (2016), as B L = A T /hw L , in which A T is the swept area of the turbine, h is the total depth of flow, and w L is the width of the local flow passage. Assuming a single cross-stream row comprising turbines of equal size and spacing, the global blockage may then be expressed as B G = n A T /hw C , in which n is the number of turbines and w C is the channel width; and the array blockage defined as B A = hw A /hw C , in which w A (= nw L ) is the array width (see Fig. 2). The local resistance is represented, following Houlsby et al. (2008) and Draper et al. (2010), by a local wake velocity coefficient α 4L , which is defined as the ratio of the velocity at the pressure equalisation point in the near wake of the turbine to the velocity measured far upstream of the turbine. In this particular code, the local blockage and resistance of the turbine(s) within a given numerical edge are determined by the average of the B L and α 4L values assigned to the two computational nodes which the edge connects (see Fig. 2). The depth and velocity of the flow passing through the turbine(s) are similarly obtained as the averages of the values calculated at these two connected nodes. The actuator disc model allows array performance to be measured using different metrics, three of which are considered in this study: the extractable power P ex , which is defined as the total amount of power removed from the flow; the available power P av , which is defined as the amount remaining when the power dissipated in local-scale mixing is subtracted from P ex ; and the global power coefficient C PG , which is defined as the ratio of P av to the kinetic energy flux of channel-scale flow measured just downstream of the channel entrance.

Validation
A non-uniform, unstructured grid is used to discretise the model domain into 17,436 triangular elements, with short sides ranging in length from 100 m in the central channel to 2 km at the ocean boundaries. A single cross-stream row of turbines with global blockage B G = 0.1 is then extended inward from one side of the channel, and a steady current is produced by establishing a fixed head difference of ζ = 0.05 m between the upstream and downstream ocean boundaries. Linear basis functions are specified and the model solutions are advanced in time using a secondorder Runge-Kutta scheme with a 1 s time step. The model is allowed to spin up from still water conditions for 2 days, after which results from the following 12 hours are extracted and time-averaged. In this steady, low roughness flow, a high turbine resistance produces a gently undulating array-scale wake, broadly similar to that which may be observed behind a bluff body at subcritical Reynolds number. The effect of this unsteadiness is to extend greatly the amount of time required for the model to achieve a (near) steady state and thus ensure that the chosen sampling period overestimates the performance of high resistance arrays. Whilst this problem is thought not to affect the overall conclusions of the present work, the effect of unsteady wakes on power production will clearly require further investigation with a more sophisticated numerical model.
Despite key differences in the underlying assumptions, the numerical array and channel models are shown to reproduce quite well the predictions from the corresponding theoreti-cal models. Estimates of the channel's maximum extractable power (also known as the channel's 'potential'), calculated for different seabed drag coefficients and using very large, full-width arrays (B L = B G = 0.8), are found to match with predictions from the theoretical channel model of Garrett and Cummins (2005) (Fig. 3a); whilst the measured variation in maximum global power coefficient with local blockage is found to agree well with results from the two-scale actuator disc theory of Nishino and Willden (2012) (Fig. 3b). The agreement with two-scale theory is shown to be poorer at high local blockage where the model is known to overesti- mate the performance of the array, but clearly sufficient to capture the leading-order physics.

Results
Optimal uniform turbine arrangements are first identified by interpolating between uniform local blockage B L and local resistance α 4L profiles. The following near-optimal uniform arrangements (to the nearest whole number of turbine edges) are then chosen as the starting points for the analysis: B L = 0.5714 and α 4L = 0.583, which produce a global power coefficient C PG of ∼ 0.966; and B L = 0.2 and α 4L = 0.424, which yield an available power P av of ∼ 2.311 MW, or ∼ 16% of the channel's ∼14.6 MW poten-tial. The arrangements which maximise P av and C PG are quite different, in this case, because the low background roughness makes the channel-scale kinetic energy flux quite sensitive to the turbine resistance (Bonar 2017). Starting with these near-optimal uniform arrangements, the effects of non-uniformity on array performance are explored by varying B L and α 4L across the width of the array. The resulting non-uniform profiles are considered bilaterally symmetric due to their reflectional symmetry about the mainland boundary from which the array is extended into the channel (see Figs. 1, 2), and include arrays with uniform B L but non-uniform α 4L , uniform α 4L but non-uniform B L , and select examples for which both B L and α 4L are non-uniform. In all cases, B L and α 4L are varied linearly across the turbine edges, with variations denoted by the extreme values at the edges nearest to (centre edge) and furthest from (end edge) the channel wall (see Fig. 2). In addition to C PG and P av , the effects of non-uniform local blockage and resistance on the local-scale extraction efficiency, which is defined as the ratio of P av to P ex , and channel-scale kinetic energy flux are also considered. For brevity, however, the corresponding figures for these additional metrics are placed in an appendix.

Non-uniform local resistance
The analysis begins by exploring the performance of arrays with uniform local blockage B L but non-uniform local resistance α 4L . Figure 4a, b illustrate the variations in global power coefficient C PG and available power P av with centre and end local resistance α 4L , normalised by the values obtained using the initial uniform arrangements.
Although the values of C PG obtained using non-uniform α 4L are not quite as high as that obtained using uniform α 4L , Fig. 4a shows that a considerable range of non-uniform local resistance profiles produce global power coefficients within ∼ 3% of this value. Figure 4b reveals a similar trend for P av , for which the corresponding range is even greater.
(The corresponding variations in local-scale extraction efficiency and channel-scale kinetic energy flux are shown in the appendix-see Figs. 10, 11.) These findings agree with those of Adcock (2015), who used a similar DG-ADCIRC model to show that although linear variations in local resistance can produce up to 5% more power in horizontally sheared flow, the performance of the array is relatively insensitive to the local resistance profile.

Non-uniform local blockage
Arrays with uniform local resistance α 4L but non-uniform local blockage B L are considered next. Figure 5a, b illustrate the variations in global power coefficient C PG and available power P av with centre and end local blockage B L , normalised by the values obtained using the initial uniform arrangements. (Figures 12 and 13 show the corresponding variations in local-scale extraction efficiency and channel-scale kinetic energy flux.) The solution space is reduced in this case because the variation in B L is subject to the additional constraint that the global blockage ratio B G must always equal 0.1. Figure 5a shows that, as with non-uniform α 4L , there is a significant range of non-uniform B L profiles which produce C PG values quite close to that of the uniform arrangement. Figure 5b reveals that the trend for P av is again similar to that for C PG and, that for both metrics, there is clear advantage to placing the larger turbines (i.e. those with higher values of B L ) on the ends of the array rather than in the centre.

Non-uniform local blockage and resistance
Lastly, four examples of arrays with both non-uniform local blockage B L and non-uniform local resistance α 4L are considered. Figures 6 and 7 illustrate the variations in normalised global power coefficient C PG and normalised available power P av with centre and end local resistance α 4L for arrays with higher local blockage B L at their ends than at their centre, and for arrays with higher B L at their centre than at their ends. (Figures 14, 15, 16, 17 show the corresponding variations in local-scale extraction efficiency and channel-scale kinetic energy flux.) The selected local blockage profiles vary from B L = 0.4 to B L = 0.7 for C PG and from B L = 0.1 to B L = 0.4 for P av . Figures 6 and 7 clearly demonstrate that the performance of the uniform arrays cannot be exceeded, or even matched, by varying both B L and α 4L together across the width. These figures also show that arrays with larger turbines at their ends consistently outperform those with larger turbines at their centre. The peak normalised C PG for the array with higher B L at its ends is not only higher than that of the array with higher B L at its centre (∼ 0.99 vs. ∼ 0.96), but is achieved with less variation in α 4L across the width (∼ 0.48 ≤ α 4L ≤∼ 0.71 vs. ∼ 0.44 ≤ α 4L ≤∼ 0.76), a higher normalised channelscale kinetic energy flux (∼ 1.14 vs. ∼ 1.01-see Fig. 15), and only slightly lower normalised local-scale extraction efficiency (∼ 1.02 vs. ∼ 1.04-see Fig. 14). The results for P av are again similar: the array with larger turbines at its ends produces a higher peak normalised P av (∼ 0.96 vs. ∼ 0.85), with less spanwise variation in α 4L (∼ 0.33 ≤ α 4L ≤∼ 0.58 vs. ∼ 0.33 ≤ α 4L ≤∼ 0.62), a higher channel-scale kinetic energy flux (∼ 0.96 vs. ∼ 0.86-see Fig. 17), and only slightly lower local-scale efficiency (∼ 1.08 vs. ∼ 1.11see Fig. 16). It is also worth noting that for all non-uniform arrays considered, power performance is maximised by tuning the smaller turbines (i.e. those with lower values of B L ) to present higher local resistances (i.e. lower values of α 4L ) than the larger turbines. Normalised C Normalised C PL TL Fig. 9 Variation in normalised: a available power P av and b applied thrust T ; across the array width W for near-optimally tuned uniform and non-uniform arrays (colour online)

(a) (b)
Normalised P Normalised T av improved by varying solely the local blockage B L or solely the local resistance α 4L across the width. For a given B L , there exists a unique α 4L to maximise either the available power P av or average global power coefficient C PG . It seems intuitive, then, that varying B L and α 4L independently of each other simply results in sub-optimal performance for the vast majority of turbines within the array. This finding agrees with Hunter et al. (2015), who used three-dimensional Reynoldsaveraged Navier-Stokes simulations of porous discs to show that the average global power coefficient C PG for a crossstream row of uniform B L is maximised by a uniform α 4L . The finding does not agree, however, with Cooke et al. (2016), who used a three-scale actuator disc model to show that a cross-stream row with uniform α 4L but non-uniform B L can produce a higher peak C PG than can be produced using the two-scale actuator disc model of Nishino and Willden (2012), for which both B L and α 4L are uniform. That being said, it should also be noted that the third scale of mixing introduced by Cooke et al. (2016) divides the single row of turbines into multiple sub-rows within the same plane, thereby creating non-uniform B L profiles much more complicated than the simple linear variations considered here.
Interestingly, it appears that neither can array performance be improved by varying both B L and α 4L together across the array width. Figures 8 and 9 show that although certain nonuniform configurations can produce similar values of C PG or P av to those of the uniform arrays, the operation of these non-uniform arrays is fundamentally different at local scale, requiring large variations in both power and thrust among the turbines, which are undesirable from a design perspective.
For all non-uniform arrays considered, the optimum strategy appears to be to tune the smaller turbines to present higher local resistances than the larger turbines. This does not appear to be a means by which to compensate for the variation in local blockage across the array width because it is not the most uniform variations in thrust and power which produce peak array performance. Rather, it implies that there is some advantage in tuning the smaller turbines, which produce the least power, sub-optimally in order to divert more flow into the path of the larger turbines which produce the most power. This would also explain why it appears to be better to place the smaller, more resistant turbines at the centre of the array rather than at its ends-because placing the smaller, more resistant turbines at the centre ensures that there are larger turbines either side to take full advantage of the flow which is diverted both left and right.

Conclusions
Simple theoretical models have provided a number of valuable insights to inform the design of tidal turbine arrays, but have only recently begun to account for the non-uniformity inherent in the flow conditions that tidal turbines experience. The work of Draper et al. (2016), which extends the classical actuator disc model to incorporate an inviscid shear flow, represents a considerable advance toward this goal. To investigate the effects of sheared flow, Draper et al. (2016) analyse the performance of an idealised turbine represented by a uniform local resistance in a non-uniform flow field. Draper et al. (2016) also note, however, that their analysis could be extended to incorporate non-uniform local resistances in order to provide a better approximation of a tidal turbine array.
In this paper, we take the first step toward such an extension by using depth-averaged numerical simulations to investigate the performance of non-uniform tidal turbine arrays in uniform flow. Results from the combined array and channel model of Bonar (2017) suggest that the performance of a tidal turbine array in an initially uniform flow field cannot be improved by varying solely the local blockage, solely the local resistance, or both local blockage and resistance together, across the array width. Certain non-uniform configurations are found to produce similar power outputs and global power coefficients to those of uniform arrays, but the operation of these non-uniform arrays is shown to require large and undesirable variations in thrust and power among the turbines. These results suggest that for an initially uni-form flow field, the optimal tidal turbine array is also uniform, that is to say that it comprises turbines of equal size, spacing, and resistance. This finding is encouraging because it is more cost-effective and much simpler to design each turbine to be the same and to operate in the same way. The result is also somewhat intuitive and, together with the findings of Adcock (2015) and Draper et al. (2016) on the performance of turbines in horizontally and vertically sheared flows, suggests a more general, and perhaps unsurprising, conclusion that tidal turbine arrays perform best when designed to account specifically for the flow conditions that they are to experience.
Finally, we note that our analysis can, of course, be extended in many different ways. The present model can be adapted, for instance, to explore the effects of varying both local blockage and resistance in non-uniform flow fields, or to incorporate channel-scale dynamics and time-variable turbine tuning strategies (e.g. Vennell and Adcock 2014;Vennell 2016). A more thorough analysis will, however, require a more sophisticated numerical model to describe more accurately the individual and collective interactions with the flow which define the performance of tidal turbines in arrays.