Effects of stator vanes on power coefficients of a zephyr vertical axis wind turbine

Renewable Energy 35 (2010) 1043–1051

Contents lists available at ScienceDirect

Renewable Energy journal homepage: www.elsevier.com/locate/renene

Effects of stator vanes on power coefficients of a zephyr vertical axis wind turbine K. Pope a, *, V. Rodrigues a, R. Doyle a, b, A. Tsopelas a, R. Gravelsins a, G.F. Naterer a, E. Tsang c a

University of Ontario Institute of Technology, Faculty of Engineering and Applied Science, 2000 Simcoe Street, Oshawa, Ontario, L1H 7K4, Canada Dalhousie University, 1360 Barrington Street, Halifax, Nova Scotia, B3J 1Z1, Canada c Zephyr Alternate Power Inc., 80 East Humber Drive, King City, Ontario, L7B 1B6, Canada b

a r t i c l e i n f o

a b s t r a c t

Article history: Received 21 February 2009 Accepted 11 October 2009 Available online 5 November 2009

In this paper, numerical and experimental studies are presented to determine the operating performance and power output from a vertical axis wind turbine (VAWT). A k-3 turbulence model is used to perform the transient simulations. The 3-D numerical predictions are based on the time averaged Spalart-Allmaras equations. A case study is performed for varying VAWT stator vane (tab) geometries of a Zephyr vertical axis wind turbine. The mean velocity is used to predict the time averaged variations of the power coefficient and power output. Power coefficients predicted by the numerical models are compared for different turbine geometries. The predictive capabilities of the numerical model are verified by past experimental data, as well as wind tunnel experiments in the current paper, to compare two particular geometric designs. The numerical results examine the turbine’s performance at constant and variable rotor velocities. The effects of stator vanes on the turbine’s power output are presented and discussed. Ó 2009 Elsevier Ltd. All rights reserved.

Keywords: Wind energy Vertical axis wind turbine

1. Introduction Due to the adverse environmental impact and non-renewable supplies of fossil fuels, worldwide efforts have turned towards renewable energy sources, particularly wind power, to meet the growing demand for electricity. Today, renewables (including nuclear power) supply about 18% of the global energy needs [1]. Through government incentives such as Ontario’s Standard Offer Program, the expansion of wind power by many private industries is growing rapidly. Although wind power dates back to ancient times, numerous technological advances over the past few decades have significantly reduced the costs of extracting energy from the wind. Although wind power has been predominantly derived from large centralized installations, there is a significant opportunity to increase capacity by smaller scale, distributed generation. Producing power where it will be used can reduce installation and land costs, with little or no transmission expense. Currently, most wind power systems consist of horizontal axis wind turbines (HAWTs), in mid to large-scale wind farms. These conventional turbines are often favored due to higher efficiency, but they are not necessarily suitable for all purposes. Horizontal axis turbines

* Corresponding author. E-mail addresses: [email protected] (K. Pope), viren.rodrigues@ mycampus.uoit.ca (V. Rodrigues), [email protected] (R. Doyle), alex.tsopelas@ mycampus.uoit.ca (A. Tsopelas), [email protected] (R. Gravelsins), greg. [email protected] (G.F. Naterer), [email protected] (E. Tsang). 0960-1481/$ – see front matter Ó 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.renene.2009.10.012

require sustained wind velocities to efficiently generate power. In an urban setting, the wind is always changing speed, and the direction is rarely uniform. In these conditions, vertical axis wind turbines (VAWTs) can be more effectively harnessed within complex urban terrains to significantly increase the capacity of small-scale wind power generation. Mertens [2] recommended three urban locations for VAWT placement. Locating it between diffuser shaped buildings can improve the aerodynamic efficiency of the turbine. Placing a wind turbine in a duct through a building can use the pressure difference between the windward and leeward sides of the building for power generation. Also, a turbine on top or alongside a building will be exposed to higher wind speed areas close to the building. Placing a wind turbine in an urban setting exposes it to significant turbulence in its operating conditions. Rohatgi et al. [3] investigated the effects of wind turbulence and atmospheric instabilities on wind turbine performance. Densely populated urban centers have unstable conditions because cities act like heat sources for the atmospheric air above. Decreased shear over the rotor diameter can be achieved by placing the turbine in low turbulence conditions. However it is generally not helpful to power generation since a stable atmospheric condition normally coincides with lower wind speeds. For complex terrain, Botta et al. [4] reported that fatigue loads on HAWTs can be increased by up to 75%, when placed in regions of high turbulence and changing wind direction. Similarly, small wind turbine loading increases by about 25% when placed in complex terrains, compared with stable wind conditions [5].

1044

K. Pope et al. / Renewable Energy 35 (2010) 1043–1051

Nomenclature A Cp D F H _ m Qevs ! r0 sr R S ! ur ! vr

vertical area of VAWT (m2) power coefficient diameter (m) force (N) height (m) mass flow rate (kg/s) equi-size skew distance to origin of rotating system (m) viscous stress (N/m2) rotor radius (m) area of mesh element (m2) whirl velocity (m/s) relative velocity (m/s)

When the change of wind direction is not considered in the overall operating conditions, the performance of wind turbines under turbulent conditions does not necessarily decrease. The effects of different levels of inlet turbulence intensity were experimentally investigated by Zhang et al. [6]. Changing the turbulence intensity from 0.9% to 5.5% or 16.2% had little effect on the turbine performance. In icing conditions, turbulence affects the heat transfer and input power needed for de-icing. Wang et al. [7] reported that the average Nusselt number for a wind turbine blade is smaller than a flat plate and cylinder at the same Reynolds number. This occurs due to differences between the turbulent boundary layer behaviour after the transition point. The authors observed that no systematic trend, relating the Nusselt and Reynolds numbers, occurred with varying liquid water content [8]. Computational fluid dynamics (CFD) is a useful design tool for wind power analysis. A large number of simulations can be performed, analyzed and optimized without investing in physical construction of many turbines with different geometrical configurations. Using CFD simulations, the torque and pressure on the rotor can be predicted. These can then be used to predict the turbine’s power coefficient. There are several CFD methods in the literature to predict wind turbine performance. One method is the Multiple Reference Frame model (MRF), which uses a time averaged solution to determine the turbine performance. It has been used frequently in wind power simulations to determine the power curves and performance of horizontal axis wind turbines [9]. Past studies have shown that numerically predicted results for a 2-D rotating rotor compare better to experimental data than predictions with a stationary rotor [10], due to the assumption of flow separation at the blade tips [10]. The Navier-Stokes equations were used with a rotating reference frame to simulate the air motion. Additional studies also found this technique gave accurate predictions of the blade forces [11]. Three-dimensional effects can lead to discrepancy with 2-D predictions [12]. Cochran et al. [13] used a 2-D sliding mesh formulation to examine Savonius turbine performance with a stator cage. A Reynolds stress 5-equation turbulence model was used and the optimal tip speed ratio (TSR) was predicted, with the use of a sliding mesh formulation to represent the rotor rotation [13]. The sliding mesh model is computationally intensive, particularly for complex geometries and it is generally limited to 2-D predictions. For a VAWT, particularly those with a stator cage, the power output can be highly erratic throughout each rotation. The sliding mesh model is time-dependent and it represents the rotor-stator interaction, thereby proving useful when predicting the turbine performance in 2-D simulations.

! v V W

absolute velocity (m/s) wind velocity (m/s) turbine width (m)

Greek

3 h q l m r s ! u U

surface roughness (m) efficiency stator angle (degrees) tip speed ratio dynamic viscosity (Ns/m2) air density (kg/m3) stator vane spacing (m) angular velocity relative to a stationary frame (rad/s) rotary speed (rad/s)

VAWTs could expand wind energy generation by utilizing an urban resource, where traditional horizontal axis wind turbines cannot operate well. Because of the attractive aesthetics of this turbine, it can be integrated into the design of buildings in a commercial, industrial or dense urban area, to produce power for onsite use. This would reduce demands on the electrical grid and environmental impact of the building design. VAWTs are omnidirectional, which negates the need for a yawing mechanism. This paper investigates the aerodynamic performance of a new type of VAWT for urban applications (see Fig. 1), called a Zephyr vertical axis wind turbine (ZVWT). This turbine has a unique design that includes stator vanes with reverse winglets (also called stator tabs). The design has a high degree of solidity, which can limit the turbine’s peak performance in strong winds. However, VAWTs with a higher degree of solidity have other performance advantages when operating at a lower TSR [14]. In particular, the effect of the tabs on turbine performance will be explored in this paper. Validation of the numerical predictions will be performed through comparisons against wind tunnel measurements. Wind tunnel experiments are presented in this paper to further examine the numerical model’s predictions of turbine performance with/ without the stator vanes. 2. Numerical formulation of air flow In this section, two fluid flow formulations will be presented: (1) a Multiple Reference Frame (MRF) model, whereby the rotor reference frame rotates with respect to the inertial frame, and (2) a moving grid transient formulation. For the second model, the overall domain is divided into two sub-domains; a stationary stator domain and a rotor sub-domain that rotates. Commercial software FLUENT 6.3.26 [15] is used for the CFD simulations. The solution domain and mesh discretizations are generated with GAMBIT 2.4.6 [16]. The governing equations, boundary conditions and solution formulations will be discussed in this section. The governing equations are given by the incompressible form of the Navier-Stokes equations. A time averaged turbulence solution is developed. This introduces several turbulent stresses into the solution. Due to the approximations needed to solve the turbulent stress equations, many different models have been developed in the past. Each method offers different advantages for the operating conditions under which it was developed [17]. The numerical predictions solve a rotating frame adaptation of the governing Navier-Stokes equations. The governing equations for fluid flow include the following conservation of mass (1) and momentum equations (2).

K. Pope et al. / Renewable Energy 35 (2010) 1043–1051

1045

pressure at the exit is estimated at standard atmospheric pressure. Standard wall functions are used to represent the stator walls, rotor blades, and outer boundary walls in the solution. The 3-D Multiple Reference Frame model is utilized for comparisons against the experimental results. To account for the high solidity and wind tunnel blockage, the numerical model is scaled to the same size (1/5 scale) as the experimental setup, and placed in a domain with boundaries at the same distance from the VAWT as the wind tunnel. The formulation uses 3 sub-domains, including a wind tunnel domain, stator domain, and the rotor domain. Above and below the rotor and stator are walls, while the other cylindrical faces in contact with the other domains are defined as interfaces. This simulates flow through the turbine. The wind tunnel domain and stator domain are kept stationary, while the reference frame of the rotor domain rotates at variable rotational speeds to represent the rotor motion. A setback with this formulation is it does not account for the rotor-stator interaction, because all calculations are based on one rotor position. The model uses an absolute velocity formulation around the interfaces between the sub-domains to provide approximate values of velocity throughout the domains. A Spalart-Allmaras model is used for turbulence effects on the MRF predictions. The Spalart-Allmaras model is a one-equation model for the eddy viscosity. The transport variable, ~v, is identical to the turbulent kinematic viscosity, except in the viscous-affected regions near the boundary walls:

2 ( ) !2 3  v~v v v 14 v  v~v 5 m þ r~v þ Cb2 r ðr~vui Þ ¼ Gv þ ðr~vÞ þ s~v vxj vt vxi vxj vxj  Yv þ S~v (3) where Gv is the turbulent viscosity production and Yv is the turbulent viscosity destruction that occurs in the near-wall region, s~v and Cb2 are constants and v is the molecular kinematic viscosity. The turbulence viscosity ratio is defined as the ratio of turbulent viscosity ðmt Þ to molecular viscosity ðmÞ as follows,

mt m

Turbulent viscosity ratioh

Fig. 1. Zephyr vertical axis wind turbine prototype – (a) illustration and (b) geometrical variables.

vr vr ¼ 0 þ V$r! vt v ! ! ! ! ! vrþ u  u  rÞ v r! v r Þ þ rð2 u  ! ðr v r Þ þ V$ðr! vt ! ¼ Vp þ Vsr þ F

(1)

(2)

The conservation of momentum contains two acceleration terms that represent the rotation. These include the Coriolis ! v r , and the centripetal acceleration acceleration, defined as 2 u  ! ! ! ! ! described by u  u  r . In these equations, r is the radial position ! from the origin of the rotating domain, u is the angular velocity of the rotor domain, ! v r is the relative velocity, p is the static pressure, s is the stress tensor and ! F refers to external body forces [15]. A mixed formulation of cylindrical and Cartesian coordinates is used to simulate the two separate regions of the domain, i.e., rotating rotor-stator section and external incoming flow, respectively. The velocity at the inlet of the domain is specified as a constant and the

(4)

The hydraulic diameter relates to the size of large eddies present in the turbulent flow. Turbulent eddies are restricted by the size of the wind tunnel, since they cannot be larger than the tunnel when the flow is fully developed [15]. The 3-D MRF numerical formulation uses the following models and parameters: finite volume method with a segregated solver, steady state, standard wall functions for near-wall treatment, wall roughness height of 0, pressure-velocity coupling by PRESTO [14]. In the next section, a comparison of CFD results with experimental data will be outlined. To allow this comparison, the majority of wind speeds simulated are obtained from the experimental test conditions. Since the experimental model uses freely moving turbine blades, as opposed to a constant rotational speed, a change in the specification of parameters is necessary. The wind speeds and rotational speeds in the CFD models correspond to those obtained in the wind tunnel. Another formulation in this paper is a 2-D transient, sliding mesh formulation. For this model, the overall domain is divided into two sub-domains (see Fig. 2 for sample mesh discretization). The rotor sub-domain rotates with respect to the stationary stator sub-domain. This formulation allows cells adjacent to the boundary between the sub-domains to slide relative to one another, thus allowing for transient predictions of the rotor interaction with the flow field, created by the stationary stator blades. To provide the appropriate values of

1046

K. Pope et al. / Renewable Energy 35 (2010) 1043–1051

C23 are C1e ¼ 1.44 and C23 ¼ 1.92. The turbulent eddy viscosity, mt, is calculated by combining k and 3 such that mt ¼ rCm k2 =3, where Cm is a constant 0.09 [15]. The model requires values of turbulence intensity and turbulence viscosity ratio. The turbulence intensity (I) is defined as the ratio of the root-mean-square of the velocity fluctuation (u) to the mean freestream velocity (uavg), described as follows,

Ih

u0 : uavg

(7)

A turbulence intensity of 10% is assumed at the velocity inlet [18] and the turbulent viscosity ratio is assumed to be 1 [19]. At the pressure outlet, a backflow turbulence intensity of 12% is assumed to account for the increased turbulence created by the turbine interaction [17]. The backflow turbulence viscosity ratio is also assumed to be 1. Other features of the solution parameters include a SIMPLE pressure-velocity coupling, second order implicit unsteady formulation, and a second order upwind scheme for turbulent kinetic energy and turbulence dissipation rate predictions. This numerical formulation will be used to predict air flow through the turbine and subsequently improve its performance. 3. Experimental wind tunnel study

Fig. 2. Sample 2-D mesh discretization of the VAWT - a) rotor, b) stator and surrounding sub-domain.

velocity for each sub-domain, a continuity of absolute velocity is enforced at the boundary. This formulation requires that the flux across the boundary between sub-domains is computed with the cell faces at the boundary misaligned. Thus, new interface zones are determined at each timestep. Each cell face is converted to an interior zone. The resulting overlap between opposing interior faces produces a single interior zone. Therefore, the number of faces at the boundary will vary with each timestep [15]. A standard k-3 model is used to predict turbulence effects in the transient predictions. This model can achieve better accuracy than the 1-equation Spalart-Allmaras model. The k-3 model is widely used and able to simulate many different flow regimes. Its two equations below are written in terms of the turbulent kinetic energy, k, and dissipation rate, 3, as follows,

v v v ðrkÞþ ðrkui Þ ¼ vt vxi vxj

"

v v v ðr3ui Þ ¼ ðr3Þ þ vt vxi vxj

#



m vk mþ t þGk þGb  r3 YM þSk sk vxj

"



(5)

#

m v3 3 mþ t þ C13 ðGk þ C33 Gb Þ s3 vxj k

 C23 r

32 k

þ S3

(6)

In this model, Gk represents the generation of turbulent kinetic energy due to mean velocity gradients, Gb is the generation of turbulent kinetic energy due to buoyancy, and Ym represents the contribution of the fluctuating dilatation to the overall dissipation rate. The variables sk and s3 are the turbulent Prandtl numbers for k and 3, with values of 1.0 and 1.3, respectively. The constants C13 and

An experimental prototype of a Zephyr vertical axis wind turbine (ZVWT) was altered to construct a 1/5th scale model with two particular design modifications, on a rapid prototype machine. The drive shaft protrudes through the anchor point and it is attached to a DC generator. This generator is connected to a variable load and a set of multi-meters for measurement. The wind speed is varied by adjusting the rotary speed of the wind tunnel fan. The measurement method is an adaptation of a previous system described by the Automatic System for Wind Turbine Testing [20]. A data acquisition method is used, with an electric motor to generate power. The system includes two multi-meters, variable resistor, and a DC electric alternator. The rotating turbine is coupled with an elastic belt to the motor to generate electricity, powering a resistor of 0.8 Ohms. Two multi-meters are connected to the circuit to measure the voltage across the resistor and current flow, thereby allowing for the calculation of the turbine power coefficient. Voltage and current readings have a resolution within 0.005 V and 0.005 mA, respectively. At any wind speed, the voltage and current vary over a small range. However, the multimeters average the values over a short time and display the results in addition to the real-time measurement. An optical tachometer measures the rotational frequency of the turbine rotor blades. The wind tunnel apparatus is a Flotek 1440 [21]. The unit has a testing area of dimensions 12 inches wide by 12 inches high by 36 inches long. The fan is located at the downstream end to generate the air flow. Control of the wind tunnel allows repeatability to within 0.25 m/s. This causes a variation of tip speed ratio between 0.02 and 0.05. The power output is determined with multi-meters to measure the voltage and current from the generator, after which power in the circuit is calculated. The TSR is measured by an optical tachometer to measure the rotary speed of the turbine. Testing was also performed to determine losses in the motor and pulley drive train. To calculate the losses in the setup, the drive train was doubled using two motors/generators, 4 pulleys and 2 belts. One motor drives the system and the second generates power. By measuring the power drawn by the first motor and the power produced by the second, the efficiency of the system was calculated. A set of batteries supplied power to the system. The power drawn by the motor was calculated by measuring the voltage across the leads of the motor and the current flowing

K. Pope et al. / Renewable Energy 35 (2010) 1043–1051

through the circuit. The motor was coupled to a 5/16th inch shaft, the same diameter as the turbine rotor shaft, by two pulleys and a belt. The shaft was then attached to a second identical motor in the same way. The power output was again measured. The input voltage was varied using two and then three batteries. This apparatus effectively doubles the drive train used in the turbine tests, therefore the losses were approximately double that of the turbine tests. Using this information, the efficiency of the motor drive train could be determined. Since the apparatus is double the actual power train, the efficiency of one half of the system is the square root of the total efficiency. This can be understood by considering each half, with the drive motor to the shaft and the shaft to the generator as separate systems, where the efficiency of both systems connected in series is the product of the efficiencies of the individual systems, ht ¼ h1 $h2 . Since each system is the same, their efficiencies should be equal and the product is h21 ¼ h22 ¼ ht . The efficiency of the drive train was found to be 45%. The power input made a significant difference in the calculated efficiency. When the input voltage is increased from 1.2 V to 2.0 V, the efficiency decreases by 3%. The output voltage and current for both supply voltages are within the range of values measured during the turbine tests, so the average efficiency is a reasonable estimate of the conditions in the turbine tests. The test setup efficiency is determined to be 38.5%. The system is more efficient at lower power levels, which suggests that the operating efficiency of each motor/generator unit differs slightly when the generator is operating at a reduced input power. There are additional sources of experimental error. The effects of scaling on the experiments are formulated by the Buckingham Pi theorem [22]. The overall problem can be represented by 8 variables, described by 3 dimensional units, so 5 pi factors are required (see Table 1). The dimensionless coefficients are the relative roughness (3/D), geometric scaling (H/D), Reynolds Number (rVD/m), TSR (UD/V), and a non-dimensional force coefficient (F/D2V2r).

Re ¼

rVD m

(8)

Matching the Reynolds number between simulations and experiments ensures that the ratio of inertial and viscous forces remain the same. Also, the 3-D MRF formulation is scaled to match the dimensions of the experimental prototype. This set of simulations was validated by reproducing the experimental conditions, such as matching the wind velocities and tunnel dimensions, to represent the wall blockage effect. 4. Results and discussion In this section, results from the numerical and experimental studies will be presented. For the 3-D mesh, the space above and below the rotor is neglected, to allow larger elements to be created. Two geometries are represented, one with the stator tabs and

1047

another with tabs removed. The full-scale turbine height (H) and diameter (W) are both 0.762 m, with a rotor radius (R) of 0.302 m. The stator spacing (s) and angle (q) are 0.266 m and 45 , respectively. The rotor and wind tunnel grids are used in all simulations, while two stator meshes are developed for comparison purposes. Two 3-D discretizations are created for each geometry: a full scale model and a 1/5th scale model. The domain is discretized with an unstructured tetrahedral mesh [23]. The quality of the mesh will be characterized by the Equi-Size Skew (Qevs), defined as

Qevs ¼

  Seq  S Seq

(9)

where S is the area of the mesh element, and Seq is the maximum area of an equilateral cell of identical circumscribing radius to that of the mesh element. A null Qevs describes an equilateral element, while a Qevs value of 1 describes an element that is completely degenerate [16]. A large number of skewed elements would result in a highly inaccurate model. The number of elements and skewness of each sub-domain is presented in Table 2. At this resolution, the results showed grid independence of the results [24]. Predictions are made using both constant and varying rotational speeds of the turbine to examine the change in performance of the turbine. One of the main practical objectives of the simulations is to determine the significance of the stator tabs. A comparison of constant and variable speeds could yield insight into the turbine’s performance with different operational characteristics. Determining the performance of the turbines at constant speed is valuable, as any power generation device connected to the electricity grid needs to operate at constant speed. The predicted optimal rotational speed is used for the simulations of the full size models. Five times the rotational speed is required for simulations of the 1/ 5th scale models. The increased rotational speed maintains a constant TSR between the full size and 1/5th scale models. Each model is operated at a variety of wind speeds to investigate how the change in velocity affects power generation. A rotational speed of 17 rad/s is maintained for the full scale models, while the rotational speed for the 1/5th scale model is 85 rad/s. To determine the power generation, the moment (M) is recorded at the bottom of the rotor where the shaft is connected and multiplied by the rotor speed (U). The power coefficient is calculated by CP ¼ M U=0:5rAV 3, where r and A are estimated to be 1.225 kg/ m3 and 0.581 m2, respectively. Comparing the coefficient of performance with TSR, the peak for the tabular design occurs at a TSR of about 0.39. With the constant operational speed, this corresponds to a wind speed of 12.5 m/s. However, the no-tab design had a peak at TSR of 0.45, corresponding to a wind speed of 10.8 m/s. Both characteristic curves can be correlated by the following quadratic curve fits (predicted by the 3-D MRF formulation): 2

Tabs : Cp ¼ 0:5285l þ 0:4422l þ 0:0059

(10)

2

No tabs : Cp ¼ 0:5674l þ 0:5225l þ 0:002 Table 1 Variables in dimensional analysis.

(11)

Table 2 3-D mesh discretizations.

Variable

Actual value

Model value

Units

Force Density Dynamic viscosity Velocity Diameter Height Surface roughness Rotary speed

– 1.225 1.82  105 12 1 1 – 1

Measured 1.225 1.82  105 Calculated 1/5 1/5 – Measured

N kg/m3 Ns/m2 m/s m m m rad/s

Sub-domain Total number of Most skewed elements element

Skewed range: Skewed range: 0.9–0.97 0.8–0.9

Tunnel Rotor Stator (tabular) Stator (without tab)

448,690 387,960 419,965

0.7708 0.9320 0.9424

0 6 27

0 33 911

606,244

0.9458

21

859

1048

K. Pope et al. / Renewable Energy 35 (2010) 1043–1051 Table 3 2-D mesh discretizations.

0.12

0.09 With tabs full size (constant load)

Cp

Sub-domain

Total number of elements

Most skewed element

Skewed range: 0–0.1

Skewed range: 0.1–0.2

Rotor Stator (tabular) Stator (without tab)

13,617 92,378

0.5971 0.5665

13,191 90,214

330 1,758

92,056

0.5469

90,278

1,508

With tabs 1/5th scale (constant load)

0.06

With tabs 1/5th scale (variable load) No tab full size (constant load) No tab 1/5th scale (constant load)

0.03

No tab 1/5th scale (variable load)

0 0

0.2

0.4

0.6

0.8

λ Fig. 3. MRF power curves with variable scaling, loading and geometry.

A peak in efficiency at a lower wind speed is useful information, since locations of lower speeds may become more suitable. When comparing the 1/5th scale model with the full size model, the curves are almost identical, which supports the notion that vertical axis wind turbines are highly scalable. When comparing the 3-D MRF tabular design with the results of the 2-D sliding mesh tabular design, the peak occurs at about the same point. The general curve is slightly lower due to the end effects that are represented by the 3-D model. Variable rotational speed operation is investigated with the 3-D 1/5th scale models. One of the limitations of the multiple reference frame model is a free spinning turbine cannot be fully simulated. Although the simulation measures the forces acting upon it, it does not react to the forces. Therefore, the rotational speed is determined by the measured rotational speed from the experimental model. With a variable rotational speed, the no-tab model outperforms the tabular model. The coefficient of performance of the turbine under these conditions is illustrated in Fig. 3. The tabular design has a predicted maximum power coefficient of 0.098 at a TSR of 0.43, which corresponds to a wind velocity of 12.9 m/s. This compares well with the constant rotational speed results. For the 2-D transient predictions, the simulation is initially conducted for 2000 timesteps with a timestep size of 0.008776 s. The

resulting 1.76 s allows four complete rotor revolutions to be completed. This ensures that a time periodic state is achieved (see Fig. 4). The moment coefficient is defined as the moment divided by 1=2rV 2 AL, where r, V, A, and L are the wind density, wind velocity, effective area of the turbine, and the length, respectively. These values were explicitly specified from the inlet conditions, with a density of 1.225 kg/m3, velocity of 12.5 m/s, area of 1 m2, and a length of 1 m. The results show a clearly defined time periodic trend. Next, the model is advanced for 500 more timesteps. This allows one more complete rotor revolution with data saved every 10 timesteps. From the results, the transient predictions of performance are obtained. For each 2-D simulation, the rotor starting position is maintained constant. The blade (identified as rotor 1 in Fig. 2) is used for all single blade comparisons. Fig. 2 shows the initial angular position of the rotor blade mesh relative to the stator sub-domain. The rotor sub-domain is removed and enlarged for clarity. For these simulations, the overall domain is discretized into triangular elements. The skewness and number of elements of each sub-domain is presented in Table 3. Grid sensitivity studies indicated that the solution converges to a power coefficient of 0.09. The selected mesh has an average discretization area of 2.191 cm2, resulting in 105,995 total cells. Comparing this grid to one with an average discretization area of 2.754 cm2 reveals a percent error of 1.9%. The results of the 2-D sliding mesh formulation show transient plots of the turbine performance. These plots are useful for determining the maximum cutoff limits, fatigue stresses and average power output of the turbine. The transient formulation provided comparable results to the 2-D MRF formulation. The predicted Cp with the sliding mesh formulation was 0.111 (a 13.3% error compared to the 2-D MRF prediction). Similar to the 3-D time averaged predictions, the

0.12 With tabs No tabs

Power (kW)

0.06

0

-0.06

-0.12 0

120

240

Angular position (degrees) Fig. 4. Time varying plot of the moment coefficient.

Fig. 5. Comparison of transient power output for a single rotor blade.

360

K. Pope et al. / Renewable Energy 35 (2010) 1043–1051

results from the 2-D transient numerical model predict a decreased performance with the tabular length design. Fig. 5 illustrates that the power output of the geometry with no stator tabs has a higher power output throughout the majority of each rotor blade revolution. Thus, the presence of the stator tabs does not offer an aerodynamic advantage and does not change the shape of the transient curves in a significant way. Although the tabular geometry exhibits the best power through the first 120 of the rotation, the following 240 of rotation has lower performance than the no-tab design.

1049

The predictions of the transient analysis provide useful results of velocity vectors and pressure contours in the flow domain. Plots generated from 50 data points were assembled to give a transient animation of the turbine and surrounding flow field. These animations can provide useful insight into the aerodynamic performance of the turbine. They can also give predictions for multiple turbine interactions. A sample of the tabular geometry is presented in Fig. 6a. From these animations, the stator tab effects on fluid flow can be visualized. The velocity vectors in the air flow

Fig. 6. VAWT a) tabular contour plot shaded by the static pressure (Pa) and b) no-tab velocity vector diagram shaded by the velocity magnitude (m/s).

1050

K. Pope et al. / Renewable Energy 35 (2010) 1043–1051

through the turbine illustrate areas of both high and low velocities (see Fig. 6b). Also, the high torque acting on the blades, as they rotate around to the backside of the turbine, can be explained by the high pressure acting on the rotor blades in this region. The area of low velocity behind the turbine is an important consideration if multiple turbines are placed in close proximity to one another. From the contour plots of static pressure, the adverse pressure on the back of the rotor blades, as they rotate into the flow stream, can be identified. This high static pressure on the backside of the rotor has a significant drain on the overall power output of the turbine. This region of the rotation explains the area of negative power on the plot of power output for one blade. Fig. 7 shows the turbine’s performance throughout one complete cycle with all five rotor blades. The plot compares the predicted transient effects for various tab angles. The base case has a tab angle of 45 . The increased tab angle geometry has a tab angle of 65 and a tab length of 1’’. The decreased tab angle geometry has a tab angle of 25 . For the three geometries, the decreased tab angle has the best performance with Cp ¼ 0.113. The base case and increased tab angle have comparable performance predictions at 0.111 and 0.109, respectively. However, the decreased tab angle geometry has the largest range of Cp values. This coincides with a sizeable range of torque forces. The increased tab angle geometry produces the smallest range in aerodynamic torque forces. This geometry would limit the degradation of electrical equipment and mechanical systems, thereby increasing the turbine longevity and power output. The predicted power output is nearly identical to the base case, while maintaining lower torque forces. If the angle of the tabs is reduced, the rotor assembly of the turbine could be subjected to destructive forces and risks of creep. During the wind tunnel experiment, data is collected for both turbine variations with 5 separate runs. A minimum constant load is used for all turbine tests. The TSR for each turbine is measured for a given wind speed. With constant loading, an increased wind speed will naturally result in an increased TSR. The curve generated by these numbers (Fig. 8) was used for a comparison of the turbine loading. The presence of loading at a given wind speed will have an effect of reducing the TSR normally observed at that wind speed, by representing a form of inertial load. With the electrical and transmission loading constant, as well as surface roughness, the variation of the curve between turbines and tests will indicate the relative amount of loading due to effects such as bearing friction, vibration

0.6 With tabs (run 1) With tabs (run 2) With tabs (run 3) With tabs (run 4)

0.4

With tabs (run 5) No tabs (run 1)

λ

No tabs (run 2) No tabs (run 3) No tabs (run 4)

0.2

No tabs (run 5)

0 0

5

10

Fig. 8. Experimental comparison of turbine performance with varying geometric design.

and shaft alignment. The variability of the curve will also indicate the repeatability of acquired data. The range of TSR shows repeatability to within 0.02 for both models. This information is further augmented with the startup and safe maximum wind speeds for the different models. The model with tabs achieved startup at a wind velocity of 8.5 m/s, and a cutout velocity of 18 m/s. The no-tab model achieved startup at a wind velocity of 7.5 m/s, and a cutout velocity of 15 m/s. A startup TSR (between 0.0 and 0.3) typically follows a linear trend. Both models displayed a tendency towards a linear decrease. For the required startup velocities, it was not possible to obtain the efficiency values for a TSR below 0.3 without the use of a breaking mechanism. The maximum measured TSR for both turbines occurred at the maximum safe wind speed, which is limited by the strength of the model. The conditions associated with higher wind speeds were mode difficult to replicate with the available laboratory equipment. The data range between 0.3 and 0.6 was typical of studies of similar turbines. Both turbines displayed a linear increase in efficiency up to a TSR of 0.5. The efficiency then peaked and began to decline again, leading to a correlated quadratic curve fit as follows (based on experiments): 2

Tabs : Cp ¼ 0:0159l þ 0:0344l 2

No Tabs : Cp ¼ 0:0126l þ 0:037l

(12) (13)

Due to a lack of data for TSR values beyond 0.6, these curves are skewed so they do not fully represent data when forecast beyond the turbine peaks. Both curves should have a steep decline after moving through the efficiency peak. For the no-tab design, the experimental results found the peak to occur at a TSR slightly higher than that predicted by the simulation, but significantly closer to that of the tabular model. The no-tab results suggest a peak of 0.12 at a TSR of 0.48, which corresponds to a wind speed of 14 m/s. The 2-D transient simulations predict a performance improvement of 13.5% with the no-tab geometry.

Cp 0.1

0.05

20

Wind velocity (m/s)

0.2

0.15

15

Base case Increased tab angle Decreased tab angle

0 0

120

240

360

Angular position (degrees) Fig. 7. Transient comparison of coefficient of performance with variable tab angle.

5. Conclusions In this paper, two independent fluid flow formulations were used to investigate the aerodynamic performance of a Zephyr

K. Pope et al. / Renewable Energy 35 (2010) 1043–1051

VAWT. The change of 2-D transient and 3-D time averaged numerical predictions of the power coefficient exhibit a comparable change in magnitude to the experimental results. This indicated that both numerical formulations provide correct trends for the changes of flow dynamics and power coefficients for changes in the VAWT geometry. The simulations indicate the scalability for both configurations. The constant speed model predicts the optimal TSR at lower wind speeds that are more likely to be found where the Zephyr turbine would be located. The simulations also exhibit a high degree of scalability for both configurations. This paper has provided useful new data and empirical correlations that characterize the operating performance of a Zephyr VAWT. Acknowledgements The authors express their grateful appreciation to the Natural Sciences and Engineering Research Council of Canada (NSERC), Zephyr Alternative Power Inc., and the Ontario Graduate Scholarship program (OGS) for their financial support. References [1] Renewable Energy Policy Network for the 21st Century (2008). 2007 Global status report shows perceptions lag reality. REN21, 75441 Paris, Cedex 9, France. [2] Mertens S. Wind energy in urban areas: concentrator effects for wind turbines close to buildings. Refocus 2002;3(No. 2):22–4. [3] Rohatgi J, Barbezier G. Wind turbulence and atmospheric stability - their effect on wind turbine output. Renewable Energy 1999;16(No. 1):908–11. [4] Botta G, Cavaliere M, Viani S, Pospo˜ S. Effects of hostile terrains on wind turbine performances and loads: the acqua spruzza experience. Journal of Wind Engineering and Industrial Aerodynamics 1998;74–76:419–31. [5] Riziotis Vasilis A, Voutsinas Spyros G. Fatigue loads on wind turbines of different control strategies operating in complex terrain. Journal of Wind Engineering and Industrial Aerodynamics 2000;85:211–40. [6] Zhang Q, Lee SW, Ligrani PM. Effects of surface roughness and turbulence intensity on the aerodynamic losses produced by the suction surface of a simulated turbine airfoil. Journal of Fluids Engineering 2004;126:257–65.

1051

[7] Wang X, Bibeau E, Naterer GF. Experimental correlation of forced convection heat transfer from a NACA airfoil. Experimental Thermal and Fluid Science 2007;31:1073–82. [8] Wang X, Naterer GF, Bibeau E. Convective droplet impact and heat transfer from a NACA airfoil. Journal of Thermophysics and Heat Transfer 2007;21(No. 3): 536–42. [9] Hahm T, Kroning J. no. 1, pp. 5–7. In the wake of a wind turbine. Fluent news, 11. Lebanon, USA: Fluent Inc.; 2002. [10] Fujisawa N. Velocity measurements and numerical calculations of flow fields in and around savonius rotors. Journal of Wind Engineering and Industrial Aerodynamics 1995;59:39–50. [11] Guerri O, Sakout A, Bouhadef K. Simulations of the fluid flow around a rotating vertical axis wind turbine. Wind Engineering 2007;31(No. 3):149–63. [12] Fujisawa N, Ishimatsu K, Kage K. A comparative study of navier-stokes calculations and experiments for the Savonius rotor. Journal of Solar Energy Engineering 1995;117(No. 4):344–6. [13] Cochran BC, Banks D, Taylor SJ. A three tiered approach for designing and evaluating performance characteristics of novel WECS. American Institute of Aeronautics and Astronautics, Inc. and the American Society of Mechanical Engineers; 2004. 42nd AIAA Aerospace Sciences Meeting and Exhibit, Reno, Nevada, Jan. 5-8, 2004. AIAA-2004-1362. [14] Paraschivoiu I. Wind turbine design with emphasis on darrieus concept. Montreal, Canada: Polytechnic International Press; 2002. [15] ANSYS Inc.. Fluent 6.3 users guide. Southpointe, 275 Technology Drive, Canonsburg, Pennsylvania 15317, USA: ANSYS Inc.; 2006. [16] ANSYS Inc.. Gambit 2.3.16 users guide. Southpointe, 275 Technology Drive, Canonsburg, Pennsylvania 15317, USA: ANSYS Inc.; 2005. [17] Gosman AD. Developments in CFD for industrial and environmental applications in wind engineering. Journal of Wind Engineering and Industrial Aerodynamics 1999;81:21–39. [18] Veers PS, Winterstein SR. Application of measured loads to wind turbine fatigue and reliability analysis. Journal of Solar Energy Engineering 1998;120(No. 4): 233–9. [19] Saxena A. Guidelines for the specification of turbulence at inflow boundaries. CFD Portal, Paris, France: ESI Group; 2007. [20] Camporeale M, Fortunato B, Marilli G. Automatic system for wind turbine testing: technical papers. Journal of Solar Energy Engineering 2001;123:333–8. [21] Flotek 1440 features and specifications. 7585 Tyler Boulevard, GDJ Inc., Mentor, Ohio 44060. [22] Munson BR, young DF, Okiishi TE. Fundamentals of fluid mechanics. 5th ed. Jefferson City, USA: John Wiley and Sons; 2006. [23] Kim SE, Boysan F. Application of CFD to environmental flows. Journal of Wind Engineering and Industrial Aerodynamics. 1999;82:145–58. [24] Pope K, Naterer G.F, Tsang E. Effects of rotor-stator geometry on vertical axis wind turbine performance. Canadian Society for Mechanical Engineers 2008 forum, June 5-8, Ottawa, Canada:2008.