Numerical study of ducted turbines in bi-directional tidal flows

ABSTRACT Installing a turbine in a duct generally increases the flow through the rotor. If there is flow confinement due to limited water depth, the presence of lateral boundaries, or the presence of adjacent turbines, a further flow enhancement through the rotor will be achieved. Extensive studies of ducted turbines in uni-directional flow have shown that rotor-based power coefficient can exceed the Betz-Joukowsky limit (BJL) for bare turbine, whereas the exceedance of the BJL is rarely reported for duct-area-based power coefficient . In this work, CFD studies have been carried out for several flanged duct turbines with different duct geometries under bi-directional flow conditions. The actuator disc approach is used to simulate the rotor. The results show that a flanged duct turbine with an inlet-arc flap and a curved brim can achieve a maximum value of close to the BJL for flow from either direction while a symmetric duct turbine- gives maximum of about 60% of the BJL. The effect of the flow confinement on shows increases with blockage ratio and a little variation in the relative performance of ducted and bare turbines. . The findings pave the way for the realization of turbines with a flanged duct.


Introduction
With the aim to reduce human-induced climate change effect and the consideration of a rapid increase in energy price, research on clean energy production is actively promoted. Global investment in renewable energy in 2017 was USD279.8 billion, a 2 percent increase over the previous year (REN21, 2018). Although wind energy is still the major source of clean energy, water flow energy (tidal energy and wave energy) has received attention, and efforts have been put into its commercialization (REN21, 2018). Turbines have been developed a long time ago and are mainly used to extract wind energy. According to Okulov and Van Kuik (2012), analytical study of turbine performance can be traced to Betz and Joukowsky in 1920. One remarkable solution of Betz-Joukowsky limit (BJL) theory is that the maximum amount of power that can be extracted from the fluid flowing through the turbine is C p = 16/27 of the total available flow power, where C p is called the power coefficient. To increase the value of C p , the concept of a ducted turbine has emerged in the past several decades. In this design the turbine is placed inside a duct such that the fluid flow will converge towards the turbine and strong vortices will be generated behind the duct resulting in a low-pressure. Lilley and Rainbird (1956) conducted a onedimensional analysis and showed that a duct could provide 65 percent more in power compared with that of a bare turbine with the same rotor diameter and the turbine performance is a function of the viscous loss, diffuser exit ratio, and the external shape of the duct. Gilbert and Foreman (1979) succeeded in increasing the fluid flow inside the duct by utilizing several flow slots to control the boundary layer but failed to prevent pressure loss by flow separation. To improve this design, Phillips et al. (1999) developed a multi-slotted duct diffuser to prevent flow separation within the duct. This device was called Vortec 7, but it failed to be commercialized due to an additional requirement of a yawing mechanism and the high cost incurred. Igra (1981) introduced the use of a ring-shaped flap to improve ducted turbine performance. His result showed that using an appropriate ring-shaped flap an 80 percent power augmentation could be obtained while about 25 percent power augmentation was achieved in its inner rear part. Hansen et al. (2000) developed and validated a computational fluid dynamics (CFD) model by modeling the rotor of a wind turbine as an actuator disk. They demonstrated that the BJL could be exceeded by a factor relative to the fluid flow rate if a mechanism was used to increase the flow rate. Lawn (2003) used one-dimensional theory to analyze the performance of a ducted turbine by treating the parts of the duct upstream and downstream of the turbine as contraction or expansion with specified diffuser efficiencies. Lawn deduced that by choosing a more lightly loaded design a ducted turbine could achieve a 33 percent enhancement in the power coefficient over a bare turbine of the same diameter. For a given duct efficiency, there is an optimum turbine resistance for generating maximum power, because the swallowing capacity of the duct increases as the resistance decreases.
Another concept of the ducted turbine was developed with a flange attached at the exit of the duct called wind-lens (Abe et al., 2005;Abe & Ohya, 2004;Ohya et al., 2008). Their computational results were in good agreement with the corresponding experimental data and showed that the wind speed inside the duct increased by 1.6-2.4 times the approaching wind speed, equivalent to about 4-5 fold power increment compared to the conventional wind turbine. Kardous et al. (2013) carried out numerical simulation, Particle-Image-Velocimetry visualization, and wind tunnel experiments to investigate the effect of the height of flange on wind velocity increase at the inlet section of an empty flanged duct. Their results showed that the presence of the flange alone caused the additional acceleration of the wind of 13 percent to 23 percent. The study, however, neglected the effect of turbine resistance in the airflow and duct interaction. El-Zahaby et al. (2017) studied the effect of the flange angle of a ducted turbine and found that flow acceleration at the entrance of the duct is optimized at an angle of 15 • to the normal of the duct with a 5 percent power coefficient increment compared to the duct with a normal flange.
The idea of ducted turbines has been extended in the context of tidal power generation probably inspired by the significant innovations in the performance of ducted horizontal axis turbines in wind turbine technology (Van Bussel, 2007). Also, Edenhofer et al. (2011) reported that horizontal axis tidal turbines are likely to follow a similar development path to wind turbines, where larger capacity turbines will be deployed depending on the progress in finding optimum duct shape (or size), blade design, and turbine array layout. The study of ducted tidal stream devices has been examined by several researchers. Setoguchi et al. (2004) investigated the performance of a turbine in a duct of symmetric shape in a wind tunnel. The duct consisted of diffuserstraight part-diffuser components and brims mounted at both ends of the duct. The case with brim height to rotor diameter ratio of 2.0 showed the best performance with an increased flow velocity of about 1.3 times the freestream flow velocity. The authors did not consider the case with a turbine present. Shives (2011) used analytical and numerical modeling to study duct of aerofoil longitudinal-section. They found that bare (non-ducted) turbine can produce more power per installed device frontal area. Shives and Crawford (2011) further conducted CFD modeling of several duct designs to gain insight into the effects of viscous loss, flow separation, and base pressure on mass flow increment provided by the duct. The authors found that viscous loss in the inlet is negligible for the cases considered while flow separation and increase in diffuser expansion, inlet contraction, and inner exit angle led to a significant reduction in turbine efficiency. They concluded that it is worthwhile to pursue further study of flanged duct turbines which create a large suction effect downstream. In placing arrays of turbines in a tidal channel, the channel confinement effect should be considered. Garrett and Cummins (2007) found that the maximum power coefficient, C p,max of a bare turbine is proportional to the blockage ratio ( ) i.e. C p,max ∝ 1/(1 − ) 2 , where = A/A c (A and A c are the rotor area and the approaching flow cross-sectional area respectively). Belloni (2013) carried out a numerical investigation of tidal turbines with symmetric-shape ducts employing the actuator disc model and CFD-BEM (Blade Element Momentum) models. She considered that C * P , the power coefficient with the outer area of the duct as the reference area, is a more appropriate parameter for comparison with the power coefficient of a bare turbine. The results showed that the performance of the devices is greatly affected by the blockage ratio and the maximum C * P for a ducted turbine is always less than the maximum C p of a bare turbine. Fleming and Willden (2016) analyzed turbines in ducts of some symmetric geometries and their results supported the work of Belloni (2013) that the maximum C * P of a ducted turbine is less than the maximum C p of a bare turbine.
For a bare turbine, the BJL caps the theoretical maximum power extracted from the fluid at about 60 percent for a single actuator disc in an unconfined environment and 0.6(1 − ) −2 for a confined environment. The BJL has been demonstrated by some researchers (Abe et al., 2005;Hansen et al., 2000;Igra, 1976Igra, , 1981 to be exceeded with the introduction of duct to generate significant suction at the duct exit and maintain a large pressure recovery within the duct (Gilbert et al., 1978). However, the demonstrations of performance gain in relation to that of a bare turbine have been questioned as the vast majority of the works used rotor area as the reference area instead of the outer area of the duct as the reference area, as well as the neglect of the effect of the blockage ratio.
The maximum values of C p (rotor area is used as the reference area) and C * P (maximum duct cross-sectional area is used as the reference area) for various ducted turbines reported in the literature are listed in Table 1. It shows that only the maximum C p s are nearly equal to or greater than the BJL while the maximum C * P s are all below the BJL. Therefore, the geometric shape of a duct for a turbine with optimum power performance remains unclear. It is reported from the literature that numerous computational, as well as experimental studies, have examined the effect of flange height, diffuser opening angle, hub ratio, and center-body length of the duct (Abe et al., 2005;Abe & Ohya, 2004;Kardous et al., 2013;Mansour & Meskinkhoda, 2014), but to the best of authors' knowledge, there has been no emphasis placed on flanged duct turbines in bi-directional tidal flows.  NaCA0015 aerofoil(deformed) = 0.0041 0.94 0.51 0.59 BJL = 0.605 Igra (1981) Shroud with ring-shaped flap (NACA 1412) n.a (pilot-plant exposed to the atmosphere) 1.79 0.45 0.74 Ohya et al. (2008) showed that the shape of the duct largely determines the streamlined velocity through the duct while Abe and Ohya (2004) pointed out that the performance of a flanged duct strongly depends on whether separation occurs within the diffuser. Thus, in this work, several carefully designed flanged duct tidal turbines with different duct shapes (symmetric and asymmetric) for both uni-directional flow and bi-directional flow will be investigated numerically to explore the power augmentation and the possibility of exceedance of BJL for C * P . Further analysis of the influence of blockage ratio on the performance of ducted turbines will also be conducted.

Modeling of the turbine
Actuator disc models have been applied to tidal current turbines (Garrett & Cummins, 2007;Houlsby et al., 2008;Whelan et al., 2009). Garrett and Cummins (2007) considered a tidal turbine in a rigid-lid channel. Whelan et al. (2009) included the deformation of the free surface while Houlsby et al. (2008) carried out a similar analysis with the inclusion of downstream mixing unlike Whelan et al. (2009). Generally, the models assume that changes in the axial velocity of the flow dominate and that crossflow velocities may be neglected. Increasing blockage improved turbine performance, resulting in higher thrust and power coefficients over a larger range of tipspeed ratios (Garrett & Cummins, 2007;Houlsby et al., 2008;Kolekar & Banerjee, 2015;Ross & Polagye, 2020;Whelan et al., 2009). Whelan et al. (2009) obtained a good approximation for a highly blocked flow pass tidal stream turbines near the free surface with the experimental tests. They concluded that agreement between the effects of free surface proximity and blockage with experiment should improve further as the blockage is reduced. Consul et al. (2013) demonstrated that the effects of rigid-lid and free-surface deformation on the power of a crossflow tidal turbine to be negligible particularly for blockage conditions less than or equal to 0.25. Kolekar and Banerjee (2015) investigated the effect of blockage on turbine performance characteristics using CFD analysis. Results showed that turbine performance generally increases with increasing the blockage ratio. No blockage corrections are necessary for blockage ratios less than or equal to 0.11 when compared to the unblocked case. Also, the performance with tip-speed ratios less than 4 was found to be weakly dependent on the blockage ratio. For instance, as the tip-speed ratio increases from 1.46 to 2.93, a 6 percent increase in C p from blockage ratio of 0.042 to 0.42 was observed. However, the impact of blockage ratio on performance curves is significant with tip-speed ratios greater than 4. The peak C p increased from 0.32 to 0.46 for the blockage ratio increased from 0.042 to 0.42. On the effect of boundary proximity, their experimental investigation suggested for optimum performance in a confined system, the turbine rotational disc should be at least one radial distance away from the solid channel wall and half radial distance below the free surface. Ross and Polagye (2020) experimentally examined the effects of blockage on the performance of a crossflow and an axial-flow turbine. Both turbines were characterized by conditions of high blockage and negligible blockage while other significant parameters were held approximately constant. Overall, the effects of increased blockage ratio increase the turbine's power performance just as in previous studies. However, error associated with blockage corrections is turbine geometries or test conditions specific.
In the present numerical simulation, the solver employed is OpenFOAM (Greenshields & Openfoam, 2016), a free Open-source CFD software package based on the finite volume method. The turbine is modeled as an actuator disc, the swirl and shedding of coherent vortices from the blades are not accounted for. This assumption may not be accurate if the turbines are close together as interaction among turbines may exist due to the swirling motion of blades. More so, the inflow crosssectional area is a circle and the actual shape of the inflow cross-section is assumed unimportant. This requires that the cross-sectional area is significantly larger than the cross-sectional area of the rotor, and the effect of the freesurface boundary should not be important. To meet these conditions a relatively small blockage ratio is required. Nishino and Willden (2012a) showed that the influence of the aspect ratio is secondary to blockage ratio when considering the efficiency of an array of turbines installed in a large rectangular tidal channel. Thus, the computation domain and the channel containing the turbine in the present simulation are assumed to be axisymmetric, which implies that the effects of free surface deformation, supporting structure and channel cross-section are negligible. The turbine region is discretized by unstructured grids around the rotor and structured grids upstream and downstream of the turbine. With the domain size normalized by the turbine diameter, D t = 20m, the inlet, and outlet were 5D t upstream and 10D t downstream of the duct respectively. Simulation results for the case with longer domain length (20D t ) downstream of the duct are less than one percent difference with the corresponding simulation results for the domain with 10D t downstream length, so the use of 10D t downstream length is deemed adequate for the present study ( Figure 1). The boundary conditions used are constant velocity, u 0 = 2.0m/s at the inlet, zero-gradient at the outlet, slip wall at the top, and no-slip wall with wall function treatment at the duct wall. For wall function treatment (appropriate for high Reynolds number flow), the dimensionless distance, y + (= u τ y/ν, where u τ is the friction velocity, y is the normal distance to the wall, ν is the kinematic viscosity of fluid) in the range 30 ≤ y + ≤ 300 was set.
The fluid is water of density ρ = 1000kg/m 3 and kinematic viscosity ν = 10 −6 m 2 /s, resulting in the Reynolds number Re = 4 × 10 7 (length scale used is the turbine diameter, D t ). The freestream turbulence intensity (FST) of 10 percent and turbulent length scale of l = 0.1D t were used in computing freestream inflow turbulence parameters. Also, a low freestream turbulence intensity of 0.1 percent was tested for the case of a bare turbine (Fleming & Willden, 2016;Nishino & Willden, 2012b). Table 2 shows the details of the simulation scheme. The turbine is treated as an actuator disc with the wake swirl and discrete blade effects neglected (Belloni, 2013;Fleming & Willden, 2016;Rahman et al., 2018;Shives & Crawford, 2011). In the numerical model, the axial resistance induced by the actuator disc was represented by a pressure drop p across the disk, where C t is a pressure coefficient or load coefficient, and where C p is the power coefficient using the rotor area as the characteristic area; and u 1 is the spatially averaged velocity at the rotor. Another power coefficient C * P can be defined using the area swept by the outer diameter of the duct A duct : Note that for bare turbine there is no duct and A duct = A, C * p = C p . The induction factor a is defined by A negative value of a refers to the condition that the velocity of flow through the rotor is accelerated and greater than the ambient velocity. In OpenFOAM, the pressure-jump p in Equation (1) is specified by using the utility function 'createBaffles'. The use of an actuator disc to simulate a turbine has the advantage that there is no need to resolve the blades explicitly and the number of grid points can be reduced significantly. Computational efficiency is important as there are a lot of scenarios to be simulated for various ducted turbines. The actuator disc imposes a streamwise resistance only and extracts the upper bound power from the flow. Additional losses due to swirl and shedding of coherent vortices from the blades are not accounted for. Furthermore, specific rotor design may be required for each duct for optimum performance. For a fair and straightforward comparison of duct performance, the assessment of the upper bound power extraction by the ideal rotor (actuator disc) is preferable. (Fleming & Willden, 2016;Shives & Crawford, 2010). All meshes were generated using Gmsh 3.0.6 (Geuzaine & Remacle, 2009) and imported into OpenFOAM using the conversion utility command, 'gmshToFoam'.
The turbulence model employed is the k − ωSST (shear-stress transport) model which is able to predict separation in adverse pressure gradient generally occurred along the entire duct surface (Menter et al., 2003;Shives & Crawford, 2011). The performance of the k − ωSST turbulence model in predicting flow separation has been extensively studied in the literature and is the most preferred turbulence model used in several ducted tidal turbine studies (Belloni, 2013;Fleming & Willden, 2016;Shives & Crawford, 2011). A comprehensive summary of the superiority of k − ωSST turbulence model in predicting complex flows was presented in Shives and Crawford (2011). The pimpleFOAM (merged PISO-SIMPLE), a robust and efficient transient solver was then used in solving the incompressible unsteady Reynolds-average Navier-Stokes (URANS) equations. A literature review on URANS is presented in section 2.3. Throughout the simulation Courant number less than unity was maintained to ensure stability. The post-processing of the results was performed using Paraview 5.3.

Choice of duct shapes
Two types of ducts are carefully chosen. The first three ducts are of asymmetric shape (Type A, B, C) and the fourth duct is of symmetric shape (Type D), as shown in Figure 2. Asymmetric shape ducts are commonly used for a ducted turbine in unidirectional flow. Their performance in bi-directional (tidal) flow is however seldom been studied. Symmetric shape ducts are always used in bi-directional (tidal) flow as their performance is independent of the direction of flow. Type B is the design used in Ohya et al. (2008) for unidirectional flow. The vertical flange at the outlet is used to increase the suction in the wake and to generate a larger flow within the duct. This duct is primarily used for unidirectional flow but its performance in bi-directional (tidal) flow will be investigated.
Type A is an improved design with a curved flange. It aims to produce better flow behavior and less head loss under reversed flow conditions. For Type C duct, a small inlet-arc flap is added to improve the inlet flow condition and reduce the inlet head loss. Type D is symmetrical and consists of a Type A duct and its mirror image. Its performance will be identical for flow from both directions. In all cases, the duct thickness, t = 0.025D t , rotor radius r = D t /2 = 10 m resulting in R/r = 1.31 (where R is the duct radius) and position of the rotor X/D t = 0. The geometric parameter of the duct is L/D t = 0.25 and the flange length to duct diameter ratio is h/D t = 0.1125. The inletarc flap length and height were set to f = 0.05D t and h f = 0.05D t , the length between the duct length and the flange periphery is S = 0.125D t (specifically for Type A and C duct types) and the inlet angle was set to θ = 4 • . The choice of these values is based on a preliminary parametric study together with the results from previous works (Abe & Ohya, 2004;Ohya et al., 2008). The blockage ratio was set at = 0.131 unless otherwise mentioned for comparisons with previous works (Fleming & Willden, 2016). Figure 3 illustrates the definition of blockage ratio ( ). In the simulation, the outer rectangular region is replaced by a circular region of the same area as the simulation is axisymmetric. In the study, the blockage ratio was kept to a maximum of 0.131 and the effect of changing the outer boundary is considered small.
The grid used in the CFD simulation is carefully chosen to satisfy the wall function requirements. A typical grid layout used in the simulation is shown in Figure 4. Grid convergence tests following the procedures recommended by the Fluid engineering Division of ASME (Celik et al., 2008) have been conducted with Type C    Table 3.
For the solutions, the order of convergence was obtained to be approximately 2. The grid convergence index for the medium grid and the fine grid was computed to be 0.95%. Figure 5 shows the sensitivity of computed C p with grid size. Since the objective is to select the most appropriate mesh that will give acceptable results in less computational time, the medium grid was used in the subsequent simulations.

URANS model with k − ωSST turbulence closure
The flow field around a flanged duct turbine is characterized as a typical bluff-body flow field with highly complex and unsteady aerodynamic and/or hydrodynamic phenomena: separation and re-attachment of the shear layer along the internal walls of the duct, vortex shedding from the flange in the external flow around the duct and mixing of the internal and external flows at the outlet of the duct. The focus of this study is to investigate the performance of flanged duct turbines under bi-directional flow and propose an optimum design. Noting that turbulence is a time-dependent 3D process whose computational simulation is particularly challenging. Among the two-equation RANS models based on Boussinesq approximation, the URANS model with k − ωSST turbulence closure was adopted in the present study. The k − ωSST model allows for direct integration through the boundary walls and shows superiority for wall simulation as previously stated. Although 3D LES (Large Eddy Simulation) or DES (Detached Eddy Simulation) is preferred to study the complex 3D flow around a bluff-body, Sun et al. (2009) demonstrated through simulation that 2D URANS solution contains 3D effect and is physical. Several other authors have extensively examined the practical applicability of the URANS k − ωSST approach. Sun et al. (2009) employed 2D URANS model with k − ωSST turbulence closure to predict the motion-induced aerodynamic forces of a rectangular section with B/D = 4 (where B is the crosssectional width and D is the section depth) at a reasonable computational cost; also the effect of the freestream turbulence was studied. A similar study was carried out by Nieto et al. (2015) using OpenFOAM. Their result also proves that 2D URANS model with k − ωSST turbulence closure can adequately model fluid-structure interaction responses, avoiding time-consuming 3D LES simulations or the development of in-house sophisticated software. Stringer et al. (2014) conducted URANS computations of flow around a circular cylinder for a wide range of Reynolds numbers. OpenFOAM was one of the solvers used. The simulation was 2D and the k − ωSST turbulence closure was employed. The computed results were close to the experimental data in terms of coefficients of lift, drag, and vortex shedding frequency at the subcritical Reynolds number range. Stringer et al. (2014) further suggested using pimpleFOAM, a solver capable of outer loop time-stepping may improve high Re convergence in OpenFOAM.
For the present problem, the flow is generated by the actuator disc which creates a pressure drop across the disc. The effect of the duct is to divert more flow towards the disc. Flow separation and vortex shedding from the duct exerts a secondary effect on the flow through the disc. The URANS model with k − ωSST turbulence closure is considered to be a viable tool for the problem. Details about the model verification are described in section 2.4.

Model verification
The first case is to simulate a bare turbine in a confined channel. The computed power coefficient, C p is compared with the corresponding power coefficient obtained from the actuator disk theory under confined flow. For this purpose, the blockage ratio of the channel is set to = 0.108 leading to C p,max = 0.7448 (Garrett & Cummins, 2007). With all flow field residuals (velocity and pressure) set to 10 −6 , the simulations were performed with a transient solver and proceeded with a time step size of 0.0001s. Results were then evaluated by comparing the time-averaged flow fields and the value of the power output coefficients (Eqs. 2 & 3). The results are shown in Figure 6. It can be observed that the computed results for the low freestream turbulence (Low FST) case match the theoretical results closely. The increase of FST to 10% (High FST) increases the strength of turbulent mixing downstream of the turbine and increases the power generation and hence the coefficient C p . In the present case with = 0.108, the increase in the maximum power output is about 8%. The present actuator disk validation results are qualitatively similar to those of Fleming and Willden (2016). The model was then applied to simulate the flow through the wind lens of Abe and Ohya  (2004), Type B duct in Figure 2, for two loading conditions: C t = 0.2 and C t = 0.9. The computed variations of pressure coefficient C Ps and normalized velocity profiles along the centerline are shown in Figure 7. The pressure coefficient is defined as C Ps = (p − p 0 )/(0.5ρu 2 0 ), where p is pressure, and p 0 is the reference pressure of the approaching flow. It can be observed that the computed centerline pressure coefficient profiles and the computed centerline velocity profiles matched well with the corresponding experimental measurements by Abe and Ohya (2004). The present k − ωSST turbulence model predicted accurately the peak pressure coefficient and the peak velocity ratio for both cases, while the k − ε turbulence model (Abe & Ohya, 2004) slightly overpredicted the peak velocity ratio and underpredicted the peak pressure coefficient for the case of C t = 0.9. In the far wake region (x/D > 3) the decreasing trend before the recovery of pressure and velocity are predicted very well by the present k − ωSST model. In the near wake region (x/D < 3) the k − ε turbulence model gives a better prediction of the profiles. The inaccuracies in the simulations may be due to the errors in the experimental measurements, the errors incurred in the turbulence models, and the errors in the grid discretization.

Power coefficients
Extensive CFD simulations have been carried out for the preliminary parametric study. The optimum geometrical parameters of the duct obtained are listed in section 2 and repeated here: L/D t = 0.25, h/D t = 0.1125, f /D t = 0.05, S/D t = 0.125, θ = 4 • . The criteria for choosing this set of parameters are based on the power output performance and the practicability in the realization of the duct. Also, the effects of duct length and flange height variations for Type B duct have been previously reported in Abe and Ohya (2004). In all cases, the rotor area is the same, and simulations are run at least for five different values of C t . The power coefficients for the ducted devices are represented by C * P and C P with different reference areas.
The power coefficient C * P is plotted as a function of the axial induction factor, a and thrust coefficient, C t for duct Type A, Type B, Type C, Type C (Reversed), Type D, and the bare devices in Figure 8. It can be seen that the asymmetric duct design (C and reversed C) show a significant improvement in power output reaching a peak value of approximately 0.8 compared to the bare device of about 0.85, whereas the symmetric duct design D attains a maximum power output of about 0.5, which is similar to that obtained from the best performing symmetrical duct designs by Fleming and Willden (2016). This significant improvement achieved by Type C duct design can be attributed to the additional flow recirculation created behind the duct device as a result of the incorporation of a flange at the exit of the duct thereby inducing a lowpressure effect (Abe & Ohya, 2004;Ohya et al., 2008;Setoguchi et al., 2004). Flow separation occurred at the tip of the flange, a low-pressure region was created downstream of the flange, which drew more flow into the duct due to larger pressure difference across the inlet and outlet of the duct. Also, the inlet-arc flap and the curvednature of the flange smoothly channel the flow into and out of the duct device. It reduces or avoids separation of flow on the near-wall region inside the duct device and downstream, in contrast to the duct equipped with a vertical flange at the exit (Abe & Ohya, 2004). The velocity field and pressure field of the flows will be further discussed in section 4.2. The optimal performance of this duct device is however achieved at lower loading conditions (C t < 1.2) and decreases with C t for C t exceeding 1.2. This is due to flow deceleration caused by the flow separation inside the duct device and the build-up of the recirculation region limiting wake expansion.
The variation of the power coefficient C p (adopting the rotor area as the reference area) with the induction factor a and thrust coefficient, C t are shown in Figure 9. The peak C p of Type C ducted turbine is about 1.6 times of that of the bare turbine while the peak C p of the Type D symmetric ducted turbine is merely close to that of the bare turbine. The magnitude of the pressure gradient increases with the size of the flange up to a certain size. Further increase in flange size will introduce a blocking effect restricting the flow into the duct (Ohya et al., 2008). Therefore, C p increases with flange size up to a certain value. Generally, the increase in velocity through the rotor is offset by the decrease in the ratio of A/A duct , resulting in a decrease of C * P with an increase in flange size.

Flow characteristics
To investigate the flow characteristics of various turbines, the computed flow streamline colored with velocity normalized by the inflow velocity at or close to the maximum power output conditions are shown in Figure 10. Similar upstream flow behavior is observed for each of the devices. Differences arise as the incident flow impinges just on the rotor and downstream of the rotor. For the bare case, the velocity decreases as the flow approaches the rotor and gradually expands further downstream of the rotor wake (Figure 10a). For Type A duct, flow recirculation from the flange occurs in the external flow. As a result, a high-speed flow develops just downstream of the rotor but is counteracted by boundary-layer flow separation along the internal walls of the duct. Further downstream, the high-speed flow quickly dies out forming a recirculation wake profile at a lower thrust coefficient (C t = 1.0) and eventually expands (Figure 10b). Similar to Type A duct, flow recirculation from the flange in the external flow occurs in Type B duct (Figure 10c). Although no separation occurs along the internal walls of Type B duct, the Type A duct demonstrates the expanding curved flange not only generates a low-pressure region behind the duct peculiar to flanged duct devices for increasing flow swallowing capacity but also align and channels the flow out of the duct, hence, a greater flow rate. This is a consequence of the increase in the aforementioned duct performance presented in Figures 8 and  9. For Type C duct equipped with a curved inlet-arc flap ensures that no separation occurs along the inner boundary walls of the duct as it streamlined the flow into the rotor plane and along the duct walls. Also, the presence of the expanding curved flange at the periphery of the duct easily align and channels the flow out of the duct while still generating a low-pressure region behind the duct. This inherent simultaneous action led to a significant increase in flow swallowing and transmission capability of the duct entailing a larger low-pressure zone (Figure 10d and Figure 11d). A higher-speed flow concentration is seen around the rotor and downstream of the duct just before the duct exit, unlike Type B with a slightly less velocity spread around the rotor. Similar behavior is observed in the reserved case (Type C Reversed) where the curved inlet-arc flap served as the duct flange (Figure 10e). This further shows that a curved flanged type is preferred and hence, the best performing duct type (see Figures 8 and 9). For the Type D duct, the front curved flange generates a large recirculation zone just outside the duct (Figure 10f). The recirculation bubble diverts more flow outside the duct and results in a lower flow rate through the rotor. It should be noted that for the ducted turbines the peak power coefficient occurred at the loading condition with the thrust coefficient significantly smaller than that of a bare turbine. Also, the incorporation of the inlet-arc flap and the shape of the inlet-arc flap and flange is a strong determinant of the nature of the flow inside and around the duct and downstream.
The pressure fields of various ducted turbines at or close to the maximum power output conditions are shown in Figure 11. For the bare case (Figure 11a), a pressure rise occurs in front of the rotor. This limits the approaching flow (at zero pressure) through the rotor. For Type A duct (Figure 11b), a slight pressure rise occurs in front of the rotor but counteracted by the action of the curved flange which diverts more flow into the duct.
For Type B duct (Figure 11c), the pressure rise is of a lesser extent due to the vertical flange which creates a lower pressure zone downstream. For Type C (Figure 11d), no pressure rise occurs in front of the rotor (unlike Type A) and a suction zone occurs just downstream of the rotor. A larger flow thus is resulted. For Type C reversed duct (Figure 11e), a slight pressure rise occurs. This is counteracted by the action of the curved flange which diverts more flow into the duct. For Type D duct (Figure 11f), the recirculation bubble just outside the duct is larger than that of Type C reversed. The downstream part of the duct with curved flange shelters the rotor and causes the low-pressure zone to occur further downstream of the rotor, resulting in a lesser flow rate.
The variation of streamwise velocity ratio and pressure coefficient along the centerline for Type C and  bare turbines are shown in Figure 12. At the location of peak power coefficient, the pressure coefficient just downstream of the rotor is the lowest, corresponding to C t = 1.05 for Type C duct, C t = 1.2 for Type C reversed duct. A peak velocity also occurs just upstream of the rotor. When the loading condition shifts away from that producing the peak coefficient, the peak velocity decreases and can be lower than the approaching velocity (case with C t = 1.2 for Type C duct), and also the value of the lowest pressure increase.

Effect of blockage ratio
In placing arrays of turbines in a tidal channel more flow will be drawn into the turbines due to flow confinement from adjacent turbines. Garrett and Cummins Figure 13. Variation of C * P with blockage ratio for bare, Type C, and Type D turbines.
(2007) have developed an analytical model to study the effect of flow confinement on the flow through a bare turbine. The effect of the blockage ratio on the peak power coefficient of the ducted turbines is studied by using CFD. Figure 13 shows the variation of the maximum power coefficient (in terms of C * P ) with the blockage ratio. As stated in section 1, the blockage ratio is defined to be the cross-sectional area of the rotor divided by the crosssectional area of the flow bounded by adjacent turbines or solid boundaries. It is observed that the maximum power coefficient increases with the blockage ratio and reaches a value of 0.78 when the blockage ratio is 0.131 (Garrett & Cummins, 2007). For a higher blockage ratio, the interaction between adjacent turbines may be significant and the inlet velocity may become nonuniform. This is beyond the scope of the study of the present work.
In the present computation, the height of the channel is varied such that the required blockage ratio is achieved. The blockage ratio varies from 0.035 (relatively unblocked) to 0.131 chosen from related past studies (Belloni, 2013;Fleming & Willden, 2016). The results for Type C and Type D turbines are shown in Figure 13. In all cases the peak C * p increases with blockage ratio. The bare turbine consistently gives the highest C * p for a given blockage ratio, as the rotor size is increased and set equal to the duct size for the other cases. The relative performance of duct turbines varies little with the blockage ratio, with Type D turbine giving the least rate. For Type C turbine, its C * p is 88% of the corresponding C * p of bare turbine over the range of blockage ratio covered. For Type D turbine, its C * p slightly varies from 62% of the corresponding C * p of bare turbine at blockage ratio 0.035 to 58% at blockage ratio 0.131.

Conclusions
A turbine with a flanged duct has been identified to give a significant improvement in power output performance for tidal flow based on CFD simulations. The duct has a curved flange at the rear end to increase the pressure drop in the wake region for forward flow and increase the approaching velocity by converging the flow into the duct. It also has a flap at the front end to reduce the inlet loss for forward flow and outlet loss for reversed flow. The maximum power coefficient C p achieved is over 60% higher than that of a bare turbine if the rotor area is used as the reference area. The maximum power coefficient C * P is close to that of a bare turbine if the largest projected area with flange is used as the reference area. The effect of blockage ratio on power coefficients is examined. In all cases, the power coefficients increase with blockage ratio while the relative performance of ducted turbines compared with the bare turbine varies little with blockage ratio.
For bare turbines, the power coefficient is capped by the BJL. Installing duct diffuser will increase C p such that the BJL is exceeded. However, C * P will generally be decreased since the effect of increasing the flow through the rotor is offset by the effect of an increase in the projected area swept by the outer diameter of the flange. The present results show that the peak C * P of a properly designed ducted turbine can approach the BJL. For practical considerations, apart from increasing power output, the installation of a duct with inlet screens can help in the protection of marine life, and the protection of the rotor.
In the present study, the rotor is simulated by an actuator disc and the power extraction is the upper bound value corresponding to an ideal rotor with an infinite number of blades. This assumption may not be accurate if the turbines are close together as interaction among turbines may exist due to the swirling motion of blades. Therefore, in the simulation, the blockage ratio is limited to the low-value range of 0.035 to 0.131. In actual application, the reduction in power extraction due to the effects of blades should be accounted for. Generally, the efficiency achieved is below 70% (e.g. Shaaban & Hafiz, 2012). A detailed study can be carried out using CFD with the blades explicitly simulated or by conducting laboratory experiments. Also, the optimum extraction of power from an array of ducted turbines in a tidal channel is a topic worthwhile for further study.

Disclosure statement
No potential conflict of interest was reported by the author(s).

Funding
The Hong Kong Polytechnic University International Postgraduate Scholarship funded by the University Grants Committee (UGC) and the Research Grants Council (RGC) of Hong Kong, SAR, China.