Simulations technique for the design of a vertical axis wind turbine device with experimental validation

Applied Energy 111 (2013) 1195–1203

Contents lists available at SciVerse ScienceDirect

Applied Energy journal homepage: www.elsevier.com/locate/apenergy

Simulations technique for the design of a vertical axis wind turbine device with experimental validation S. Rolland ⇑, W. Newton, A.J. Williams, T.N. Croft, D.T. Gethin, M. Cross College of Engineering, Swansea University, Swansea SA2 8PP, UK

h i g h l i g h t s  Numerical techniques for the efficient transient simulation of vertical axis wind turbines.  Comparison with experimental data with excellent agreement under baseline conditions.  Limitations of the simplified model highlighted, in particular regarding turbulence modelling.

a r t i c l e

i n f o

Article history: Received 16 April 2012 Received in revised form 23 March 2013 Accepted 8 April 2013 Available online 29 June 2013 Keywords: Computational fluid dynamics Vertical axis wind turbine Experimental validation

a b s t r a c t Computational fluid dynamics (CFDs) simulation tools are developed to analyse the aerodynamic performance of a novel design vertical axis wind turbine (VAWT) and compared against careful data from a 1.6 m diameter device in a wind tunnel. The study investigates the extent to which a CFD model employing the simplest turbulence representation can provide a useful input to evaluate the impact of several key operational parameters: wind speed, rotor speed, yaw angle and blade pitch angle. The results show that simple turbulence modelling techniques are sufficient to evaluate the performance of the turbine in the designed operating conditions and can predict when the turbine will run outside each parameter’s operational range. However, such a simple turbulence representation may need further refinement to quantify more precisely the extent to which performance is affected when running some way outside this range. Ó 2013 Elsevier Ltd. All rights reserved.

1. Introduction Wind energy is a growing industrial sector fuelled by economic, environmental and political motivations [1]. The latest global trend report published on the Global Wind Energy Council website showing figures updated to 2010 describes how in the previous decade, the global installed wind capacity grew twelve fold to reach a production capacity of 197 GW with nations such as China doubling their capacity annually [2]. The economic impact of this growth, under the impulse of global policies, means that the wind turbine industry is now a key sector for employment and investment. Research work is ongoing in both academic and industrial centres to improve the efficiency and reliability of the wind turbine devices that will service the growing energy demand. In the present work, numerical simulation techniques, based on Computational Fluid Dynamics (CFDs) software technologies, are used to investigate the response of a vertical axis wind turbine to design changes with the additional benefit of experimental validation. The same approach was taken by Chong et al. to calibrate compare ⇑ Corresponding author. Tel.: +44 01792606871. E-mail address: [email protected] (S. Rolland). 0306-2619/$ - see front matter Ó 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.apenergy.2013.04.026

numerical results and experiments in H-rotor designs, for example [3]. The device considered here is a simplified version of a vertical axis wind turbine (VAWT) designed in collaboration with an industrial partner, with a view to be relevant to the industry while simplified enough to be published without compromising the company’s intellectual property. Prior work has also been shown to be of value to compare numerical and experimental design with early stage prototypes yielding low efficiency [3]. Previous work has shown that VAWTs can work under conditions which are not normally favourable for traditional horizontal axis wind turbines (HAWTs) [4], thus generating interest from some sections of the energy generation industry sector [5,6]. The work presented below makes use of CFD based simulation technology to model experimental benchmark tests on the novel design VAWT as illustrated in Fig. 1 [7]. A finite volume code, PHYSICA, was used to perform the CFD simulations. The ultimate aim of the work is to provide a practical simulation tool, which can be of use to the wind turbine industry, thus offering the possibility of rapidly examining a large number of design configurations. Due to the nature of the VAWT considered here a number of challenges were encountered. First, the turbine was not axi-symmetrical. This implies that transient simulations were required to reproduce the

1196

S. Rolland et al. / Applied Energy 111 (2013) 1195–1203

Fig. 1. Illustration of the studied VAWT design, perspective and 2D projection.

variations of the turbine output, unlike HAWTs, which can be simulated effectively using symbolic representation in steady state simulations, see, for example, the work of Williams et al. [8]. To solve this challenge, a sliding mesh technique was implemented in the CFD code. This enables a geometrically exact simulation of the rotor’s motion with time. As with many other simulation techniques for rotating machinery, the origins of the sliding mesh technique lies in simulations of mixing in stirred chemical reactors [9– 12] and aero engines [13,14]. The device studied here operates on an impulse or drag principle, rather than lift [1] and no publication relating numerical and experimental work was found in the literature, whereas lift driven devices have recently been the focus of this type of study [15]. The work here describes the technique used to simulate the wind turbine and the assumptions used to build a simplified CFD model. The results from the simulations are then examined quantitatively and qualitatively against the measurements and observations from the set of wind tunnel experiments by Rolland et al. [7]. The experimental device is a very simple partially faired vertical axis wind turbine. The diameter of the turbine is 1.59 m. The blades’ chord is 0.25 m. The fairing covers 90° inside and outside the turbine. The inner part of the fairing is almond shaped defined by to circular arcs of radius 0.54 m. The outer fairing has a circular inner face with a 0.8 m radius and extends out, such that the outer face is a tangent to the rotor 2 m in length. 2. Computational modelling strategy The strategy to capture the behaviour of the turbine involves the mathematical modelling of the wind load on the rotating structure. Thus, the physical behaviour to be captured involves the movement of a viscous fluid with the properties of air. The governing conservation equations, which represent the transport of momentum and mass of a Newtonian fluid, when expressed in vector notation are:

@ ðquÞ þ r:ðquuÞ ¼ r:ðlruÞ  rp þ S @t

ð1Þ

@q þ r:ðquÞ ¼ Sm @t

ð2Þ

Here q, u, l and p represent the density, velocity, viscosity and pressure of the fluid. S represents a source of momentum and Sm is a source of mass. In the scenario considered here the mass source, Sm, is zero. The flow rates are slow enough that the fluid can be considered as incompressible so that the transient term, @@tq ; in Eq. (2) is also zero. Hence, this mass conservation equation simplifies to:

r:ðquÞ ¼ 0

ð3Þ

The flow regime around the turbine is expected to be turbulent. This is normally reflected through an additional mathematical model to capture a measure of the turbulence – these effects are then reflected through their impact upon the viscosity of the fluid [16]. There is a challenge in dealing with rotating flows, and most practical one- and two-equation turbulence models do not capture the impact of the non-isotropic nature of the rotating flows, and so become of very little value in a context suitable for the present case. Moreover, the flow is naturally transient and the most effective way to embed the effect of turbulence was to utilise an enhanced viscosity over that of laminar. That is, the turbulence capture could be considered as a simplified LES model where the viscosity is simply specified rather than calculated through some sub-grid model [16]. If the mesh is fine enough then the transient nature of the calculations should enable the larger eddies to be resolved well enough to make reasonable predictions of the torque and power output generated under a range of wind loading conditions. This approach has proved to be quite effective in a number of other scenarios involving swirling or rotating flows [17]. To predict the fluid flow around the turbine, it is necessary to solve these equations for the unknown variables, velocity and

Fig. 2. Illustration of the conforming sliding mesh technique.

S. Rolland et al. / Applied Energy 111 (2013) 1195–1203

1197

Fig. 3. Mesh structure with boundary condition faces highlighted in red.

pressure fields, with respect to initial and boundary values within a space and time domain. A numerical solution is obtained through CFD simulation technology. The technique employed here to model the turbine is the cell-centred finite volume method extended to an unstructured mesh in three dimensions [18–20]. The technique relies on a discretisation of the simulated domain into small control volumes or cells. The software used in the present case is Physica [18], which uses the standard finite volume cell centred discretisation scheme referenced above. Additionally, the RhieChow interpolation scheme [21] is applied to the momentum equation (1) to minimise spurious oscillations of velocity in the computational solutions and the SIMPLEC method [22,23] is used to ensure that the pressure field is coupled with the continuity equation (3). One of the challenges specific to the simulation of the turbine is the transient simulation of the rotor motion. Unlike some pumps or stirred reactors, for example, the geometry surrounding the rotor is not axisymmetrical, thus rendering impossible the use of techniques with computationally implied motion such as the rotating reference frame technique; see for example, Croft et al. 2009 [24]. Instead, the motion was modelled explicitly throughout the time stepping scheme with the use of a conforming sliding mesh technique. The discretisation grid is generated with two sub-domains, one is fixed while the other one rotates and they share a cylindrical boundary. Each face of the elements on the shared boundary exactly maps to another face on the other side of the boundary, as illustrated in Fig. 2 below. In the context of this work, a single mesh is created and then split on the interface between the static and moving parts in order to create duplicate faces and points on the interface. In order to reduce the problems associated with elements changing volume the mesh divisions on the split are set to give the same area for all faces on the interface. All this is done in such a way as to enable the CFD code to run with scalable performance on parallel clusters [25]. To measure turbine performance, the torque on the blades was computed and it can be used to calculate the power output of the turbine. The power output was computed as the sum of torque contributions from the faces on the surface of the blade. These torque contributions were obtained as the cross product of the force on the face, based on pressure p, area dA, centre x, and normal to the face n, with the radius vector to the axis of the rotor, r. This tor-

que vector, s, is projected onto the rotor axis, defined by the vector D through a point xo, to preserve only the effective torque component, a scalar value dT.

3. Simulations and results 3.1. Power output The simulations were designed to validate the above CFD model against experimental data of a simplified VAWT system. They were therefore set-up to reflect the wind tunnel environment in which experimental data was acquired [7]. Assumptions were made regarding the oncoming flow to the device. The first one is that the device is high enough to render the end effects at the top and bottom plates of the turbine as negligible. This assumption was tested and validated in the experimental work published concurrently [7]. This led to the simplifying assumption that the flow could be simulated using a two-dimensional model – that is, a slice through the plane normal to the axis of the VAWT system, as illustrated in Fig. 1. It was also assumed that the oncoming flow is uniform at 20 rotor diameters upstream of the simulated turbine and the pressure returns to its baseline value 105 diameters downstream of the turbine. The simulation domain width is 7.94 m, as is the width of the test section of the wind tunnel used to obtain the experimental data. Non-slip boundary conditions are applied at the domain walls [26]. An inlet is set up where the air velocity boundary condition is applied and the outlet is specified as a pressure boundary condition with a reference pressure p = 0. Finally, it was assumed that a direct resolution of the larger turbulence scales and the use of an enhanced fixed viscosity was sufficient to take into account turbulent losses. This relies on the fact that a restricted range of wind velocities is applied to the turbine and the driving mechanism is based on impulse, thereby using the drag properties of the blades rather than lift

Fig. 4. Mesh split in eight partitions for parallel computing.

1198

S. Rolland et al. / Applied Energy 111 (2013) 1195–1203

Fig. 7. Simulated and measured power output up to 20 m/s.

Fig. 5. Results from the study on the impact of artificially enhanced kinematic viscosity.

Fig. 8. Simulated and measured power output with varying blade pitch angle.

Fig. 6. Simulated and measured power output against rotor speeds under 8 m/s and 11.5 m/s.

was mapped as a function of the rotor speed expressed as the tip speed ratio (TSR, k). The TSR is the expression of the rotor speed with respect to the wind speed calculated as:

k ¼ Xr=U properties, which is less demanding to capture in a CFD model. Consequently, a relatively fine mesh was used and an enhanced air viscosity value was used. The key question then is – what is a suitable value for the viscosity? The spatial discretisation was carried out by splitting the mesh into two sub-domains to use the conforming sliding mesh technique. The static sub-domain is composed of two regions: the inside of the rotor, which contains static elements of housing, and the outside of the rotor. The rotor itself is contained in the rotating sub-domain, and annulus. This led to the coexistence of two sliding interfaces on the inside and the outside of the annular rotating sub-domain, as illustrated in Fig. 4, which shows the core geometry and the way in which the mesh is partitioned. A typical simulation involved a mixed mesh with some 140,000 elements; 52% of the elements are located within the rotor (Fig. 3). Three full rotations of the rotor are made during a simulation in 1440 time steps, by which time a stable flow regime had been established. All the results are then based on the data from the final rotation. The code was run on a parallel cluster at Swansea University, which comprises 24 nodes each having 2 INTEL Nehalem quad core processors yielding a potential total of 192 processes, which have a 2.93 GHz clock rate. The interconnect is a Myrinet system with a latency of 110 ns and a full duplex bandwidth of 10 Gb/s. The software is run under Linux CentOS 5.5 and the communication is managed by the MPICH library. A typical simulation was run on 8 processors and took 43 CPU-hours. Once a mesh independent simulation condition had been determined some tests were run to see if a single viscosity could be identified to best match the experimental data. The power output

ð4Þ

in which X is the rotational speed of the rotor in radians per second, r is the radius of the rotor and U is the wind speed. The tip speed ratio is used as a common measure of rotor speed independent of wind speed [1]. In Fig. 5 the experimental simulation results for the power output are shown for an oncoming wind velocity of 8 m/s in the default configuration [2], at 30° to the wind and with a blade pitch angle of 15°. It is clear that for a viscosity of 0.004 m/s2 (i.e. 200 times the laminar value), an excellent match was obtained with the measured values of power output. This single viscosity value was used in all the other simulations reported here. It is recognised that while this is not ideal from the perspective of a very accurate capture of the turbulence phenomena, the study shows that it adequately reflects the performance of the device with respect to the power output. It is clear, from the extended discussions in the literature, exemplified by that in Davidson [16], for example, the practical capture of the impact of turbulence is not straightforward, especially for complex flow scenarios, such as those considered here. One reasonably successful route has been identified by Tucker [27–29], for flows with lift structures, who used an LES model in the body of the flow field and a two-equation model in the neighbourhood of the internal structure to capture the details of the boundary layer. The work here concerns a rotating drag device where the boundary layer flow is of no real significance and so the approach taken was to try what might be regarded as a zero-order LES model – that is, the flow field was solved the with no explicit turbulence model and assumed a simple enhanced value

S. Rolland et al. / Applied Energy 111 (2013) 1195–1203

1199

Fig. 9. Layout of the pressure tappings around the turbine base.

Fig. 10. Comparison of the average of experimental pressure values and the simulated average.

for the viscosity. Fig. 5 shows the sensitivity of the simulated power output to variations in the viscosity value as compared to the measured data. Clearly, more energy needs to be dissipated as power output is over-estimated by the model when low viscosity values are used. There is little sensitivity to the value at low viscosity levels as the sub-grid turbulence dissipation is not fully accounted for. However, as the viscosity level is increased, the power output sensitivity increases. Based on the blade length, the Reynolds numbers in simulations associated with the different values are ranging from 105 down to 102. The value of viscosity found to match experimental

data yields a Reynolds number of 5.102, therefore returning an essentially turbulent flow to a regime, which can be solved as a laminar flow. This choice involves monitoring the performance of the model throughout the modelling exercise against the experiments and its limitations are highlighted below. Four parameters were studied to compare numerical and experimental results based on the power output of the device: the sensitivity to wind speed, rotating speed of the rotor, angle of wind incidence, named yaw angle in this study due to the adoption of the wind tunnel terminology, and the blade pitch angle.

1200

S. Rolland et al. / Applied Energy 111 (2013) 1195–1203

Fig. 11. Comparison of the experimental (left) and simulated (right) pressure curves for all tested wind speed at tip speed ratios of 0.3.

Fig. 12. Experimental and simulated pressure curves with blade angles of -25°, -15°, -5° and 5°.

Fig. 13. Simulated and measured power output at yaw angles between 0° and 50°.

The experimental values were measured at the MIRA full-scale aerodynamic facility. All values were calculated and corrected as laid out in the companion publication to the present work [7] and were based on the guidelines from MIRA [30]. The first set of tests considers the power output of the device under two specific wind speeds: 8 m/s and 11.5 m/s. The test

results are plotted in Fig. 6 at wind speeds of 8 m/s and 11.5 m/s. The comparison of numerical and experimental results demonstrates excellent correlation. The simulations were run at a range of tip speed ratios and wind speeds to evaluate performance under a broad range of wind loadings. The predicted and measured turbine power outputs are shown in Fig. 7 with good correlation. Throughout the explored range of wind speed, the optimal tip speed ratio for power output was found to be 0.4. The excellent correlation between the experimental and simulated power output under a range of operating conditions does not translate to all variations in the turbine configuration. Fig. 13 shows the power output from the turbine at varying yaw angles, which is equivalent to changing wind direction. As explained above, the numerical model is set up for direct resolution of turbulence. With calibration, it is possible to reproduce the behaviour of the wind turbine, however, turbulent phenomena arises that cannot be resolved, in particular within the chosen resolution of time and space discretisation. The choices, based on practicality, have the effect of altering the predicted response of the turbine compared the measured response. The discrepancy in Fig. 13 between predicted and measured results is interpreted as a limitation of the turbulence modelling strategy. Although

S. Rolland et al. / Applied Energy 111 (2013) 1195–1203

1201

Fig. 14. Comparison of the experimental (left) and simulated (right) pressure curves for all tested wind incidence angles at tip speed ratios of 0.3.

Fig. 15. Turbulent vortex pattern developing between the blades of the turbines comparison between experiment and simulation.

Fig. 16. Turbulent vortex pattern developing between the blades of the turbines comparison between experiment and simulation.

the experimental evidence in [7] does not suggest a drop in performance, the pressure profile from the simulation shows that there is a significant change in the flow when the wind incidence angle is 40° or greater. This is further reviewed in the section on model limitations. The impact of the blade pitch angle on the power output could not be simulated very accurately. The experimental results, Fig. 8, show that the noise in the data is comparable to the power response of the turbine to the change in the blade pitch angle, with the numerical results showing that the peak performance is

predicted at 25° but measured at 15°. The experimental data points show that the repeatability at 15° and 5° is within 12%, which extended to the points with no repeated experimental results yields simulation results within experimental accuracy. 3.2. Pressure measurements In the experiments, pressure tappings were used at a number of locations around the turbine. The locations of sensors on the experimental turbine are documented in the published experimental

1202

S. Rolland et al. / Applied Energy 111 (2013) 1195–1203

Fig. 17. Comparison of the flow patterns between the experiments and the simulations as the vortex is shed from between the blades into the main wind stream.

results [7] and shown in Fig. 9 for clarity. The results from the simulation were obtained by retrieving coordinate-based pressure values from post-processing software. Both experimental and simulated pressure values are variations from a baseline pressure value set to zero in the simulations and set as atmospheric pressure as measured at the time of each experiment in the case of the experimental values. Because the pressure reading times in the experimental setup cannot be realistically claimed as instantaneous (2 s reading time), the pressure values from the simulations are calculated as the average pressure over a revolution of the rotor for comparability. Additionally, the experimental results are shown as the average reading over several experiments. Data on experimental repeatability is available from [7]. The pressure profiles show that there is an excellent correlation between the experimental and the numerical results, showing that the flow patterns are adequately captured by the CFD model. High and low pressure regions are found where expected. A more detailed analysis shows that the maximum and minimum pressure values are slightly overestimated by the model, in particular at sensor references F5, L1, M1 and M10 (Fig. 10). The readings at sensors L1, M1 and M10 suggest that the CFD model predicts that the high pressure region is on the back face of the housing, whereas it is measured to be closer to the tip, maybe even on the face inclined towards the rotor. The pressure tapping F5 is located diametrically opposed to the housing. It is in the location where vortex shedding initiation was noted. Fig. 11 shows that the differences between computed and measured pressure values appear to scale with the wind speed. However, although the reference pressure was taken into account, no clear relationship could be established to correlate these differences. It is clear from Fig. 10 and Fig. 11 that there remains an element of discrepancy between the baseline pressure in experiments and simulations. It could be taken into account by shifting the curves so that they show the same average pressure. This was attempted and no relationship could be found. In addition, Fig. 11 shows that the correlation between experimental pressure readings at different heights of the turbine (series L and M) do not present any significant differences. Other figures are therefore simplified to include only the M series of pressure readings from the experiments. The response to the change in blade pitch angle is not measurable as a change in pressure field, as illustrated in Fig. 12. Additional sensors in different locations may have been necessary. Fig. 8 corroborates the fact that the blade design results in little sensitivity to its pitch angle within the range studied here as the simulations also show little sensitivity to the blade pitch angle. 3.3. Limitations of the current modelling technique Although the model is performing well in the chosen working conditions, limitations have been found. A study was performed

to measure the response of the turbine to a change in yaw angle. This type of study may be necessary on a design to evaluate the impact of gusts from varying directions, for example, or to ensure the integrity of the turbine in case of breakdown of the yaw control mechanisms. It was found that the response to a changing yaw angle is not accurately predicted by the model. Fig. 14 illustrates the main weakness in the economical simulation strategy selected: the pressure values are predicted to span a larger range than measured. Additionally, the response to varying yaw angles is perceptibly altered. The numerical model predicts significant pressure building up on the rear face of the housing from a yaw angle of 20°, whereas the measurement show that this does not happen until 30°, with the 10° increment resolution. Additionally, the pressure on the upstream corner of the housing (measurement M10) is predicted to peak at 30°, whereas it is shown to increase throughout the range experimentally. The pressures along the edge guiding air to the rotor are also showing dramatic reversal of the pressure patterns from 40°, indicated that the numerical model predicts complete flow separation. A considerable pressure drop was indeed measured but to a much lesser extent. 3.4. Smoke tests A smoke wand was used in the wind tunnel to visualise the flow patterns. The recording equipment to capture the results consisted of four video cameras mounted on the base of the turbine. The images below show the turbulent flow patterns developing on the open side of the turbine, diagonally opposed to the housing in the area above sensor F5, Fig. 9. Because their position and direction were clearly established, it is straightforward to locate everything observed on film; moreover, the cameras were placed in such a way that the perimeter close to the turbine was captured. The camera situated at the forward corner away from the housing yielded very good images. It shows the development of turbulence between blades. The three dimensional aspects, such as counterrotating features (Fig. 15 and Fig. 16) do not appear on the twodimensional simulations, of course, but the growth of the turbulent features with the rotation as well as its release into the wind stream is correctly predicted (Fig. 17). 4. Conclusions The ultimate objective of this work is to develop a CFD based computational model which provides a useful tool for the assessment of novel candidate designs for VAWT devices. This paper summarises the development of the computational model, which is based upon a conventional unstructured mesh finite volume CFD code and uses as simplistic approach to capturing the effects

S. Rolland et al. / Applied Energy 111 (2013) 1195–1203

of turbulence. This approach is a crude LES type model where the viscosity is fixed a some enhanced level (with respect to its laminar value) by tuning the simulation to a subset of the experimental data and then using it to validate CFD model against the remaining data. The experiments covered a wide range of operational parameters and the CFD model was carefully compared with data for: a) The power output. b) The pressure distribution around the device. c) Videos of smoke illustrating features of the flow such as vortex formation, etc. The CFD model reflects the experimental data well overall and numerical values are shown to give very good correlation within the working range of each studied parameter with respect to overall power output and also the pressure distribution. Notable deviations from the experimental results were observed in the CFD model when extreme values of geometrical angle were used. The authors attribute these deviations to the simplified turbulence modelling approach. Two phenomena are non-physical in the chosen approach: for some configurations at the outer end of the operating parameter range, an excessive amount of energy is dissipated through fluid shear and the flow is artificially maintained in a laminar regime under conditions, which would otherwise be turbulent. The authors acknowledge a need for further work on turbulence modelling for this class of problems. Having identified this need, the next phase will investigate the discretisation requirements for a direct simulation and the use of two-equation, etc. models and compare with a view to identifying the best solution. Acknowledgements The authors wish to acknowledge the Welsh Government and C-fec Ltd. for the funding and support throughout the Project, notably Messrs Anthony Fenwick-Wilson and Kevin Mowbray for their input in this Project. References [1] Hau E. Wind turbines – fundamentals, technologies, application, economics, 2nd ed., Berlin Heidelberg: Spring-Verlag; 2006. [2] Global Wind Energy Council. Global wind report. Belgium: Brussels; 2010. [3] Chong WT, Fazlizan A, Poh SC, Pan KC, Hew WP, Hsiao FB. The design, simulation and testing of an urban vertical axis wind turbine with the omnidirection-guide-vane. Applied Energy; in press. [4] Eriksson S, Bernhoff H, Leijon M. Evaluation of different turbine concepts for wind power. Renew Sustain Energy Rev 2008;12:1419–34. [5] Balduzzi F, Bianchini A, Carnevale EA, Ferrari L, Magnani S. Feasibility analysis of a Darrieus vertical-axis wind turbine installation in the rooftop of a building. Appl Energy 2012;97:921–9. [6] Chong WT, Naghavi MS, Poh SC, Mahlia TM, Pan KC. Techno-economic analysis of a wind–solar hybrid renewable energy system with rainwater collection feature for urban high-rise application. Appl Energy 2011;88(11):4067–77.

1203

[7] Rolland S, Thatcher M, Newton W, Williams AJ, Croft TN, Gethin DT, et al. Benchmark experiments for simulations of a vertical axis wind turbine. Applied Energy 2013;111:1183–94. [8] Williams AJ, Croft TN, Masters I, Bennet CR, Willis MR. A combined BEM-CFD model for tidal stream turbines. In: Third international conference on ocean energy, Bilbao, Spain; 2010. [9] Pericleous KA, Patel MK. The source–sink approach in the modelling of stirred reactors. Physico-Chem. Hydrodyn. 1987;9:279–97. [10] Bujalski W, Jaworski Z, Nienow AW. A CFD study of the homogenization with dual Rushton turbines – comparison with experimental results. Part II: the multiple reference frame. Trans Inst Mech Eng – Part A 2002;80:97–104. [11] Jaworski Z, Bujalski W, Otomo O, Nienow AW. CFD study of homogenization with dual Rushton turbines – comparison with experimental results. Part I: Initial studies. Trans Inst Chem Eng – Part A 2002;78:327–33. [12] Ng K, Fentiman NJ, Lee KC, Yanneskis M. Assessment of sliding mesh CFD predictions and LDA measurements of the flow in a tank stirred by a Rushton impeller. Trans Inst Chem Eng – Part A 1998;76:737–47. [13] Sayma AI, Vahdati M, Imregun M. An integrated non-linear approach for turbomachniery forced response prediction. Part I: formulation. J Fluids Struct 2000;14:87–101. [14] Vahdati M, Sayma AI, Imregun M. An integrated non-linear approach for turbomachnery forced response prediction. Part II: Case studies. J Fluids Struct 2000;14:013–125. [15] Howell R, Qin N, Edwards J, Durrani N. Wind tunnel and numerical study of a small vertical axis wind turbine. Renew Energy 2010;35:412–22. [16] Davidson PA. Trubulence: an introduction for scientists and engineers. Oxford, UK: Oxford University Press; 2004. [17] Humphreys NJ, McBride D, Shevchenko DM, Croft TN, Withey P, Green NR, Cross M. Modelling centrifugal casting: the challenges and validaiton. In: Proceedings of the international symposium in liquid metals processing and casting, Nancy, France; 2011. [18] Croft TN, Pericleous K, Cross M. PHYSICA: a multiphysics environment for complex flow processes. Numer Methods Lamin Turb Flow 1995;3:1269–80. [19] Chow P, Cross M, Pericleous K, Koulis A. A natural extension of the conventional finite volume method into polygonal unstructured meshes for CFD application. Appl Math Model 1995;20(2):170–83. [20] Croft TN. Unstructured mesh finite volume algorithms fr swirling turbulent reacting flows (PhD thesis), London: University of Greenwich; 1998. [21] Rhie CM, Chow WL. Numerical study of turbulent flow past and airfoil with trailing edge separation. AIAA J. 1983;21:1525–32. [22] Patankar SV, Spalding DB. A calculation procedure for heat, mass and momentum transfer in three-dimensional parabolic flows. Int J Heat Mass Transf 1972;15(10):1787–806. [23] Van Doormal JP, Raithby GD. Enhancement of the SIMPLE method for predicting incompressible fluid flow. Numer Heat Transfer 1984;7:147–63. [24] Croft TN, McBride D, Cross M, Gebhardt JE. Multi-component free surface flows and rotating devices in the context of mineral processing. Int J Comput Fluid Dynam 2009;23:93–107. [25] Croft TN, Carswell D, Cross M, McBride D, Rolland S, Slone AK. Parallel computational fluid dynamics – not without its challenges. In: Proceedings of the 1st international conference on parallel, distributed and grid copmuting for engineering, Pécs, Hungary; 2009. [26] Versteeg HK, Malalasekera W. An introduction to computational fluid dynamics – the finite volume method. 2nd ed. Harlow, UK: Pearson Eductation Ltd.; 2007. [27] Tucker PG. Novel MILES computations for jet flows and noise. Int J Heat Fluid Flow 2004;25(4):625–35. [28] Tucker PG. Computation of unsteady turbomachinery flows: Part 1-Progress and challenges. Prog Aerosp Sci 2011;47:522–45. [29] Tucker PG. Computation of unsteady turbomachinery flows: Part 2: LES and hybrids. Prog Aerosp Sci 2011;47:546–69. [30] MIRA Ltd. Registered in England No. 402570. MIRA aerodynamic wind tunnel facilities – users’ handbook. Nuneaton, UK; 2007.