Numerical investigation of the wake interaction between two model wind turbines with span-wise offset

Wake interaction between two model scale wind turbines with span-wise offset is investigated numerically using Large Eddy Simulation (LES) and the results are validated against the experimental data. An actuator line technique is used for modeling the rotor. The investigated setup refers to a series of experimental measurements of two model scale turbines conducted by NTNU in low speed wind tunnel in which the two wind turbines are aligned with a span-wise offset resulting in half wake interaction. Two levels of free-stream turbulence are tested, the minimum undisturbed level of about Ti ≈ 0.23% and a high level of about Ti ≈ 10% using a passive upstream grid. The results show that the rotor characteristics for both rotors are well captured numerically even if the downstream rotor operates into stall regimes. There are however some difficulties in correct prediction of the thrust level. The interacting wake development is captured in great details in terms of wake deficit and streamwise turbulence kinetic energy. The present work is done in connection with Blind test 3 workshops organized jointly by NOWITECH and NORCOWE.


Introduction
Nowadays, the two constraints of space and economics prevents the location of wind turbines sufficiently apart from each other.Therefore, the inevitable wake interactions between wind turbines occur which means that the downstream wind turbines operate in the wake of the upstream ones.The wake of a wind turbine is characterised by stream-wise velocity deficit and high turbulence intensity.The wake deficit results in less available power for the downstream turbines and as a consequence of high turbulence level, the structural loads increase which shorten the life span of the rotor blades.It is reported that losses in energy production due to wake interference effects are 10 − 20% depending on the layout of the wind farm [1].The detailed unsteady wake development and wake interactions are studied by researchers, and their common goal is cost reduction through optimizing the efficiency of wind farms and minimizing the loads.
The prediction of the wind turbines performance as well as the wake development are approached by various means such as experimental techniques and numerical models.The experimental investigations are conducted rarely on large (full) scale and often on small scale wind turbines.Acquiring reliable experimental data from full scale wind turbines is usually difficult and expensive because a complete test set needs a full documentations of the incoming free-stream velocity and turbulence.However, limited well-documented full scale investigations are available, such as the NREL Phase VI rotor measurement at NASA Ames wind tunnel [2] and the MEXICO project using the DNW wind tunnel [3].The provided data has served as a reference test cases for several numerical predictions such as [2; 4; 5].Due to lack of space in the full scale measurements the detailed study of far wake development as well as wake interactions are conducted using small scale model wind turbines [6][7][8].The experiments are conducted using variety of small scale rotor sizes, different airfoils used along the rotor span, operating conditions and layouts.These studies can provide useful information about some main features of the wake structure.There is however, disadvantages due to the low Reynolds numbers and other scaling effects.
Recently, Adaramola et al. [9] performed a series of experiments on 0.9 m diameter model wind turbines.First, they measured the performance and near wake profiles of a single model wind turbine [10].Wake interference effect was then studied on the performance of the downstream wind turbine, with two turbines aligned.The work was further extended by considering the half wake interaction condition with the downstream turbine shifted span-wise [11].
This report presents numerical simulations of the wake interaction between two model wind turbines with span-wise offset using a technique that combines an actuator line method with large eddy simulations (LES).Simulations are carried out for two free-stream turbulence levels.The predicted wind turbines performance and their wake development are then compared with the experimental data [11].The objective of the work is to contribute to an improved understanding of interacting wakes and to study how the external aerodynamic blade loads of a wind turbine is affected upon impingement of a wake of an upstream turbine.

Numerical modeling
The Navier-Stokes computations are carried out using the 3D flow solver EllipSys3D [12; 13].The code solves the discretized incompressible Navier-Stokes equations using a finite volume approach.In the computations, the flow field around the wind turbine rotor is simulated using the actuator line technique [14], where the presence of rotor in the computation is replaced by the body forces.The forces along the blade are computed using the local flow fields and tabulated airfoil data.
The computational domain is a regular Cartesian grid, divided into a number of blocks.The domain is composed of approximately 24.5 million mesh points with each blade represented by 43 actuator line points.High resolution grid points are needed in order to capture gradients and to resolve the wake.The present grid is chosen from a study performed in [15].In the computations, the time step is chosen such that the tip of the blade cannot pass more than one computational grid points at each time step.Therefore, the time step is set to ∆t = 0.003.The applied boundary conditions are constant inflow velocity, convective outflow and the wind tunnel walls are included using a slip Euler boundary conditions.The wind tunnel turbulence is modeled by introducing synthetic resolved turbulent fluctuations generated with the Mann algorithm [16].
The computations are carried out using model wind turbines used in NTNU experiments, see [17].The modeled rotors are three-bladed horizontal axis with 14% thick NREL S826 airfoil along the span, operating at an effective sectional Reynolds number around Re=10 5 .In the computations, the rotor performances are directly related to the lift and drag performance of a rotor blade and therefore reliable blade data plays an important role in the result accuracy [18].A measurement campaign of the S826 airfoil was therefore undertaken in order to get reliable airfoil characteristics at the relevant Reynolds numbers.The measurement of the 2D airfoil performance and the experimental condition are reported in details in [15].The results show that the profile is quite sensitive to Reynolds number variation and the blades are subjected to static stall hysteresis.The higher lift and lower drag values are obtained at increasing incidences whereas, the lower lift and higher drag coefficients are measured at decreasing incidences.The higher level of lift and drag coefficients are used for the rotor computations and they are reported in figure 1.To resemble the experiments, the wind turbine towers are also included in the numerical simulations.The towers are modelled in the same manner as fixed actuator lines in which the body forces are used instead of the real geometry of the towers.Here, the local forces along the tower (F cyl ) consist of lift and drag forces where C L and C D are the lift and drag coefficients, ρ is the air viscosity, u is the local velocity, d is the tower diameter.e L and e D denote the unit vectors in the direction of lift and drag forces, respectively.Schewe [19] measured mean drag and rms lift coefficients of a smooth circular cylinder for different Reynolds numbers.Based on his measurements, the drag coefficient is considered as a constant C D = 1.2 while the lift coefficient has a harmonic variation with time computed as C L = A sin(ωt) where A = 0.3 and ω = 0.4π is the wake shedding frequency.t is the non-dimensional simulation time.

Numerical setup
The numerical investigations are conducted on two modelled wind turbine, which are aligned with a streamwise separation of ∆x = 6R and span-wise shift of ∆y = 0.83R.The wind tunnel layout is shown in figure 2. The rotors are positioned in such a way that the projected area of the upstream rotor covers half of the downstream turbine.There are small differences in the rotor diameters and tower shapes of the two modelled wind turbines.We will refer to the report by Krogstad et al. [11] for more details about the experiment.Two levels of turbulence intensity are investigated, minimum level around (T i = 0.23%) which corresponds to measurements of the empty NTNU wind tunnel and a high level of turbulence with T i = 10% using a passive upstream grid to consider the effects of atmospheric turbulence.The synthetic resolved turbulent fluctuations is introduced in a cross-sectional plane at approximately x/R = 1.In the measurements, the free-stream turbulence intensity levels are measured at the position of the first rotor when the rotor is removed.A similar approach is followed for the computations, in which we first introduce the turbulence fluctuation to the domain in absence of the actuator lines and then we compute the turbulence quantities.Figure 3(a) shows that the predicted free-stream turbulence levels are in a good agreement with the reported wind tunnel conditions.As there is no significant span-wise variation in the flow, it is expected that the turbulence intensity decays slowly in stream-wise direction.This decay rate shown in figure 3(b) implies that generally a good agreement is achieved.

Power and Thrust Characteristics
The performance of the downstream turbine is investigated with the upstream turbine operating at its design tip speed ratio (λ where Ω 1 is the angular velocity of the upstream turbine.The performance of the downstream rotor is computed over a range of tip speed ratios λ 2 .The power and thrust coefficients of the rotors are computed as where P and T are the power and thrust forces, respectively.Note that the coefficients for the downstream turbine is estimated using the uniform inlet velocity of U ref = 10 m/s.The performance envelope of both rotors in presence of the low free-stream turbulence is reported in figure 4(a-b).Considering the upstream rotor (full symbols) we can observe that the power output increases with the tip speed ratio, reaches a maximum value at the design tip speed ratio and then decreases with further increase in tip speed ratio.The first indication of the stall is observed at λ 1 = 4 and with further decrease in λ 1 the blade operates in deep stall over the entire span.At very high tip speed ratios of λ > 11, the near root region operates at negative angle of attack which results in negative power production.The model wind turbine shows a rapid drop in power production for λ 1 < 4.This is most likely related to the static stall hystereses of the blades.The effect of the interacting wake on the downstream turbine is also presented.At the first glance, a general reduction in production level of the downstream turbine is observed.This is due to the fact that the downstream turbine is exposed to a relatively lower velocity respect to the upstream one.The numerical prediction is in excellent agreement with the experiments for the power coefficients.The thrust coefficients in figure 4(b) shows a steady increasing rate by increase of the tip speed ratios.The thrust coefficients for the upstream and downstream turbines are approximately the same, implying that a similar physical force is applied to the rotors, however, the available energy for the downstream rotor is reduced.Here we can observe that the thrust coefficient of the upstream rotor is predicted well by the numerics, while the thrust level of the downstream rotor is under-predicted at all the tip speed ratios.
The performance of the rotors are also investigated in the presence of the high free-stream turbulence.By comparing figure 4(c) with figure 4(a) we can observe that the increase in free-stream turbulence results in a decrease in the power production of the first rotor at the optimum rotor performance and the hysteresis phenomena is also removed in the experiments.Once again, a good agreement is achieved between numerical simulations and the experiment.Figure 4(d) shows that the thrust coefficients of the second turbine increases by almost 50% for λ < 4 and reaches 30% for the runaway condition.The thrust coefficients of both upstream and downstream rotors are not in agreement with the experiments as there is a 23% over-prediction of the thrust for the upstream rotor and considerable under-prediction for the downstream one.

Wake Properties
To investigate the wake development of the interacting wakes, three conditions for the downstream rotor are examined, the partly stall condition (λ 2 = 3.5), optimum power performance (λ 2 = 4.75) and run-away condition (λ 2 = 8.0).The former condition occurs where the downstream rotor operates into its stall regime.In the latter condition, the rotor blades are partly operating at a propeller state.
The interacting wake development at optimum rotor condition is shown in figure 5.At low free-stream turbulence, the wake behind the first turbine is seen to be dominated by organized tip and root vortices.The interaction between the root vortices and the tower wake results in the root vortices break up immediately after the tower.Afterwards, this unsteady wake starts to interact with the tip vortices resulting in tip vortex breakdown.The flow velocity is reduced significantly as it passes by the second rotor.The deceleration of the flow is partly a result of the rotors being enforced to operate at pre-defined tip speed ratio.The generated wake behind the second rotor is asymmetric and it is dominated by small-scale turbulence structures.Further downstream, we can observe that the wake center tends to move toward the centre of the wind tunnel.Figure 6 shows the wake development by introducing high free-stream turbulence.Here, in addition to very unstable root vortices the tip vortices are also broken down fast due to the instability mechanisms triggered by the presence of the high free-stream turbulence.The same mechanisms are also observed for the two other rotor conditions (λ = 3.5, 8.0) with slight differences.Two quantities inside the wake behind the second rotor are computed to quantify the structure of the flow field, the time-averaged stream-wise velocities and the stream-wise turbulence intensity.Then, two horizontal profiles are extracted for each quantity at specific streamwise positions one close to the rotor e.g.x = 12R and the other one further downstream at x = 16R.The wake deficit profiles for different rotor conditions are shown in figure 7. By looking at the profiles for low free-stream turbulence, we can observe that both experiments and predictions produce asymmetric wake.Comparing the wake deficit profiles at x = 12R, it can be observed that regardless of the operating condition of the downstream rotors, the wake size and its shape are correctly predicted by the numerical simulations, although some differences in the wake profiles are observed.In the experiments, the rotors are placed inside wind tunnel at span-wise positions of y 1 /R = −0.45 for the upstream rotor and the downstream rotor at y 2 /R = +0.45.Here, we can observe that the profiles in i.e. figure 7  regions, the region (−1.5<y/r<− 0.5) where the wake is only affected from the upstream rotor, the central part of the wake (−0.5<y/r<1.0)which is a blend of the two wakes, and, lastly the region (1.0<y/r<1.5)where only the effect of extracting the energy from the flow by the downstream wind turbine is observed.In run-away condition (figure 7 (e)), the outboard of the rotor blade has a high local thrust while the inward part operates at propeller state.This results in relatively high deficit near the tip of the blade and relatively low deficit near the root.Small differences in the wake deficit of the near wake can be related to the difference in prediction of the thrust coefficients; in addition the presence of hub is not accounted in the simulations.In the far wake, the diffusion modifies the flows significantly so we observe Gaussian wake deficits, indicating that the flow here is dominated by small-scale turbulence structures and the wake shows very small details of its source.
The effect of turbulence on the near wake deficit is almost negligible while the main differences are observed in the far wakes.There is generally excellent agreement inside the wake area for all the considered conditions however, a general disagreement exists between the computations and the experiments outside of the wake region, where the speed up occurs due to presence of the wind tunnel walls.This can be explained by the fact that in computations the boundary layer development over the walls are not considered.
The normalized stream-wise turbulent stress profiles are shown in figure 8.The turbulent kinetic energy distribution in the near wake is characterized by four distinct peak generated from the tip vortices.It is observed that these peaks are smoothed out by the span-wise turbulent diffusion further downstream.The numerical predictions are in agreement with the experiments and the position of the tip vortices are correctly predicted.

Conclusions
The performance of the interacting wakes between two rotors with span-wise offset has been investigated numerically.The CFD computations are conducted using large eddy simulations (LES) and actuator line method.The rotor computations are conducted using two free stream turbulence levels.The numerical results are then validated against the wind tunnel measurements.
The presented results show that the power losses for a turbine operating in the wake of another turbine is significant.In addition, the increase in free-stream turbulence results in the decrease of the power production of the downstream wind turbine.The general shape of the power and thrust coefficients of the downstream rotor are correctly predicted but the thrust coefficient is consistently lower than the experiments and the thrust coefficient for the upstream turbine is over-predicted.The under-prediction of the thrust coefficient increases by increase of the free-stream turbulence.
Analyzing the wake structure close to the downstream rotor shows that the wake deficit profile is asymmetric and the profile consists of three regions, upstream rotor wake, downstream rotor wake and the overlapping wake region.Further downstream, the diffusion modifies the flows significantly and this three regions are smoothened and we can observe the Gaussian shape deficit.The prediction of the wake deficit as well as the turbulence intensity profiles generated from the half wake interaction of the two model wind turbines are in reasonable agreement with experiments and the differences is mostly caused by the difference in thrust coefficients.
On the basis of the information provided, we have demonstrated that complicated LES computations can be a complementary tools for investigation of the wind turbine wake interactions.

Figure 1 .
Figure 1.a) Lift and b) drag coefficient distributions of NREL S826 at different angle of attacks.The reported Reynolds numbers are in 10 3 .

Figure 2 .
Figure 2. Schematic view of the wind tunnel layout.The inlet plane is at x/R = 0 and the flow direction corresponds to the positive x-axis.The span-wise distance between the rotors are 0.83R.The horizontal lines, A-A and B-B, show the position of the extracted profiles.All the parameters in the simulations are normalized using the blade radius R and the inlet velocity U ref .

Figure 3 .
Figure 3. a) Computed time averaged turbulence intensity profiles obtained along horizontal line through the rotor center.The average low level background turbulence is around 0.23%.b) The turbulence intensity decay profile taken along the wind tunnel for high background turbulence condition.The experimental measurements at different stream-wise positions are shown in open circles and the numerical prediction in the line.

Figure 4 .
Figure 4. a) Power and b) thrust coefficients of the two modelled wind turbines in presence of low free-stream turbulence.c) Power and d) thrust coefficients of the two modelled wind turbines in presence of high free-stream turbulence.The measurements are presented as the black circles and the computations as the red triangles.The filled symbols indicates the upstream model turbine and the open symbols for the downstream one.The uncertainty in power and thrust coefficients of the experiments are estimated to be less than ±4% and ±2%, correspondingly.

Figure 5 .Figure 6 .
Figure 5. Contours of a) vorticity magnitude, b) averaged stream-wise turbulence intensity, c) instantaneous axial velocity and d) time averaged streams-wise velocity in a x/y-plane cut show the development of interacting wake for the optimum condition (λ 2 = 4.75) and low freestream turbulence.The flow reaches its fully developed conditions after non-dimensional time of t = t.(Uref /R) = 65 equivalent to approximately three flow through domain and the solutions are time-averaged over non-dimensional simulation time t = t.(Uref /R) = 42 equivalent to two flow through the domain.

Figure 7 . 2 refFigure 8 .
Figure 7. Mean velocity profiles U /U ref computed along horizontal line through the rotor center.(a,c,e) represent the conditions where low free-stream turbulence is considered and (b,d,f) are for the high free-stream turbulence.(a-b) λ 2 = 3.5, (c-d) λ 2 = 4.75 and (e-f) λ 2 = 8.0.The uncertainty in measurements of the mean velocity is estimated to be less than ±4.6% U ref .