On the wind speed reduction in the center of large clusters of wind turbines

Journal of Wind Engineering and Industrial Aerodynamics, 39 ( 1992 ) 251-265

251

Elsevier Science Publishers B.V., Amsterdam - - Printed in The Netherlands

On the wind speed reduction in the center of large clusters of wind turbines Sten Frandsen Wind Engineering Section, Department of Meteorology and Wind Energy, Riso National Laboratory, DK-4000 Roskilde, Denmark

Summary The paper intends to probe the feasibility of a method to determine the wind speed reduction in the center of a large wind farm by means of simple boundary layer theory. Therefore, the simplest possible assumptions have been chosen. The paper presents an approach for determination of the reduction of the wind speed assuming the wind turbines effectively act as roughness elements. The model makes use of similarities to so-called canopy flows, where the surface drag and the drag on individual obstacles are added to form the total drag. Results are compared with existing models for reduction of efficiencywithin wind turbine clusters and good agreement is found.

1. Introduction

Consider the wind speed at hub height in a large wind farm far away from boundary effects, i.e. so far from the boundaries of the wind farm that hub height wind speed is constant. Under these circumstances, there must be horizontal homogeneity when disregarding the local disturbances of the machines, so to make up for the extraction of energy by the wind turbines and for turbulent dissipation, a vertical energy transport is the only source of fresh energy. The basic idea is that the flow characteristics mean wind speed and root mean square of fluctuations in wind speed within the wind farm may be separated in a horizontally averaged component and a component dependent on the local conditions (i.e. near-wake conditions):

u(X) fUm(X)+u~(x) and a(X)fam(X)+aw(X),

(1)

where u and a are mean wind speed and standard deviation on wind speed at the position X, respectively. Thus, the spatial averages of Uw and aw are zero. In this paper focus is on the spatial average values Umand am, which deep inside the wind turbine cluster are assumed to be functions only of height above ground. Before presenting the model the basic concepts on an early wind farm model

0167-6105/92/$05.00

© 1992 Elsevier Science Publishers B.V. All rights reserved.

252

developed at Riso and some elements of atmospheric boundary layer theory are given. 2. Background:

momentum

theory

The model to he presented is not dependent on the validity of the simple Betz theory. However, for the testing of the model it is advantageous to maintain Betz theory because of the simple relations between the relevant parameters. Therefore, some main features of the Betz theory is summarized here.

2.1. Betz theory Characteristics of the cylindrical Betz stream tube are summarized in Table 1. Note that for convenience the flow speed downstream is set to ul - Uo( 1 - a) and the flow speed through the rotor to ul - Uo(1 - ½a), where a is the induction factor. Employing Bernoulli's equation and the continuity equation for the stream tube alone give the rotor thrust on the "idear'rotor and power from the wind turbine:

T=½pu2Ar2a(1-½a) and P=½pu3oAr2a(1-½a)2,

(2)

where p is air density, Uo free-stream, velocity, and Ar swept rotor area. The quantities CT=2a(1--½a)

and

Cp=2a(1-½a) 2

(3)

TABLE 1 Characteristics of Betz's stream tube. The parameter a is the axial induction factor: the flow speed immediately downstream is u l = Uo(1 - a ) , where Uo is the upstream flow speed Parameter

Upstream

Rotor disc

Downstream

Flow speed Rotor area

uo Ao

u, = Uo( 1 - ~a) A , - 1 -Ao-~a

u i = Uo( 1 - a) Ao A l - 1-Za

Rotor area

Ao

Ar

Rotor diameter

ro

Rotor diameter

ro

r,=~ rr

1 - ½a

ro

A1--A, i _ a ro r,--~½ a rx= r,~] 1--a

253

are the thrust and power coefficients 1, respectively. The power coefficient has maximum for a = 2/3: max (Cp) = 16/27.

2.2. Cluster efficiency: momentum theory A number of models for the cluster efficiency of wind farms have been devised, e.g. refs. [1-3]. Here, the model developed by Jensen [4] will be used for comparison. The model is simply based on the assumption that momentum is conserved in the wake, see Fig. 1, which leads to the following flow speed in the distance x from the wind turbine in the wake downstream:

u-'u°[1--a(l+o~-(x/rl))]' 1

2

(4)

where Uo is the upstream flow speed, a is the induction factor in flow speed immediately behind the wind turbine and c~ is the momentum entrainment constant. Note that the initial radius of the wake, rl is the downstream radius of the Betz' stream tube, related to the rotor radius by rl - r~x/(1 - ½a) / (1 - a). Employing the superposition principle on the wakes of an infinite row of wind turbines, Fig. 2, each with the same flow induction a = ~, the asymptotic flow speed reduction factor can be shown [4] to be

R,_U_~_ l Uo

2y 1-y'

a(

1

Y=2\l+a-(x/rl))

~2

"

(5)

By more or lessempirical means, the entertainment constant has been estab-

r-- o¢X + r l rl

Fig. 1. Schematic wake model with definition of symbols. tin conventional Betz formulation the flow speed down stream is expressed as u~ = Uo(1- 2a'). This makes the expression for Cp=4a' (1 -a' )2 with m a x i m u m at a' = I/3.

254

Fig. 2. Schematic multiple wake model.

lished as a = 0.5/ln (h/zo) '. Here, h is hub height of the wind turbines and Zo is the surface roughness of the area of the wind turbine cluster. 3. Background:

ABL

Some basic concepts of the Atmospheric Boundary Layer theory to be employed herein are summarized in the following. 3.1. Basic concepts, no external shear forces

-

3.1. I. Neutrally stratified wind profile Neglecting the displacement height, the vertical mean wind profile is found in neutrally stratified atmosphere - to be

u . - / ~ In

,

(6)

where u. is the friction velocity, k is the yon Kfirmfin universal constant (equal to 0.4), z the height considered, and zo the surface roughness. Experimental determination of each of the constants in Eq. (6) turns out to be difficult. The expression is usually assumed valid from z ~ 10Zoto the height of a few hundred meters. 3.1.2. Shear stress For a boundary layer in balance the effective shear stress is independent of height: r.=pu~ .

(7)

255

The stress in the flow is equal to the flow friction with the terrain surface. Likewise, the standard deviation of fluctuations of the along-wind component, u, is not a function of height: a. ~ 2.5u..

(8)

3.1.3. Vertical energy transport and dissipation In a horizontally homogeneous shear layer, there can be only vertical net transport of energy. The vertical power flux through I m 2 is r(z)~(z), ref [7]. The input power at the top of a column with unity cross section is equal to the dissipation in the column plus the output power at the surface: z

z

Pinm~

~(Z')(lz'

~t-Pout==~ ~

~(Z t )

o

dz'-ru-r.u.--r.[u(z)-u.].

(9)

o

By differentiation, the dissipation per

m 3

d u. 1 pu3. e ( z ) = r . - ~ ( u ( z ) - u . ) = r . k z - kz "

is found: (10)

Assuming that Reynold's number (defined as Re =La./v, where L is the scale of turbulence and O is the air viscosity ( = 1 . 5 X 10 -5 m 2 s - l a t 2 0 ° C ) is constant during the first period of decay, the dissipation may also be expressed [71 as

(11)

e=PAaL ,

where A is a constant. This expression is all equivalent to the expression of Eq. (10). Assuming the macro scale to be l ~ 20z [8] equalizing Eqs. (10) and (11) and introducing the value of au from Eq. (8) give an evaluation of the constant A: A ~ 3.2.

3.1.4. The geostrophic drag law The wind speed at high altitude is linked to the friction velocity by the geostrophic drag law: u.

G

G = In(.~,,~o) ~Zo = ~ exp

,

(12)

where f is the Coriolis parameter ( = 1.2 X 10 -4 at latitude 55 ° ), and G the geostrophic wind speed. The geostrophic drag law may be used to eliminate the roughness parameter from the log-profile equation: u.

-~ln

.

(13)

256

3.2. Shear flow with obstacles For evaluation of the flow characteristics in between the roughness elements of the atmospheric boundary layer the commonly accepted tool is to assume that the forces from the individual obstacles may be evenly horizontally distributed.

3.2.1. The mean profile Consider a vertical column of air with unity cross section. Employing horizontal projection we find the forces on the column are shear forces on the top and the bottom, and a distributed load, t (z), caused by drag from some obstacles. Thus, the load from the obstacles acting of 1 m of the air column is

(14)

t = --pct u2 ,

where ct is the drag coefficient per m height, per m 2 land occupied by the obstacles:

c,(z)-l~., ~Di(z)ci(z) --At

[m-']

(15)

where Di (z) and ci (z) are the cross-wind dimension and drag coefficient of obstacle No. i, respectively, and At is the total area occupied be the obstacles. Horizontal balance of forces leads to the following equation: Z

r(z) rc-f t ( z ' ) - -

dz' =U=~'~=pctu ^ dr

2,

(16)

c

where z - c is the constant height of the bottom of the column. In the so-called mixing length formulation 2 - frequently used in the research of the flow within crop canopies, e.g. refs. [9-11] - Eq. (16) form the basis for a numerical solution to the flow problem in between the roughness elements, in our case the wind turbines. 2Mixing Length Formulation:A popular assumption concerning the horizontalshear force isthat it may be expressed as r =pl du/dz),where I isthe so-calledmixing length. Insertingthis expression for r in Eq. (19) gives an equation to li.nkthe mixing length and flow speed: d(ldu/dz) 2

~ ,2[du~ d2u

. dl ~du~.~

=½c, u

By certain assumptions about the mixing length, this differentialequation may by means of numerical methods provide a solutionto the flow problem in question.

257

3.2.2. Vertical energy transport, turbulence

Consider again a vertical air column with unity cross section. The energy inflow a the top less the energy out-flow at the bottom of the column will equalize the integrated dissipation plus possible energy extraction: Z

Z

S

r(z)u(z)=rcUc=~ ~(z') dz' + pCpU3(z')dz" ..

("

(17)

¢

where ~ is dissipation and Cp is the distributed power extraction coefficient defined similarly to ct. Differentiating Eq. (17) with respect to z yields dr du u -~ -b r - ~ = ~-b pCpU3 .

(18)

Inserting Eq. (16) in (18) yields 3

du

(19)

~--'p(ct--Cp)U -br dz"

When the first term on the right side is zero one arrives at the usual expression for dissipation in the free stream: ~ - r du/dz. Obviously, this would imply, that ct-cp-O.

Equalizing the two expressions for dissipation, Eq. (11), in Eq. (19) yields Ap--~--p(ct--Cp)U a~_ . - 3+ r-~=~ du L ,/s

(20)

3 r

Knowing L as functions of height, and the vertical wind profile makes it in principle possible to calculate the vertical turbulence profile. 4. Model applied

The weakness of the mixing length model referred to previously is that it need as input a suggestion for the variation of the mixing length (and thus the scale of turbulence) with height. For the very high density of obstacles in canopy models [ 9 - l l ] that way of modelling leads to applicable results, assuming that the mixing length is constant: the wind profile is found to be exponential within the crop canopy. In the case of a cluster of wind turbines the density of obstacles is typically much less and the height of the obstacles much larger than for the c r6p canopy. Also, experiments [ 12 ] indicate that below the "rotor layer" - the layer swept by the wind turbines rotors - the vertical wind speed profile is reduced but still logarithmic.

258 4.1. Mean flow Thus, in the case of a wind farm one cannot expect the mixing length to be constant. To avoid making a priori decision on the scale of turbulence the following crude assumptions are made: (1) The impact on the flow of the wind turbines may be considered evenly, horizontally distributed. (2) The change in shear due to the impact of the wind turbines is concentrated at hub height, Fig. 3. (3) The spatially averaged vertical wind profile is logarithmic over and under hub height. (4) The shear is likewise constant over and under the "rotor layer" (z < h - r , or z > h + r,) different by t, (5) It is assumed that the wind turbines are evenly spaced the distance x, so that the dimensionless separation is s - X / r r . Thus, one wind turbine occupies the area (Srr) 2. With these assumptions the change in shear stress per m 2 land may be expressed as h2

t--r.2--r.2= f petu2 dz ,

(22)

hz

%2

ro~

-r

--T,2

"--,

--,oCtUh2

--%1 . / / Wind Speed Fig. 3. Illustrationof verticalflow shear and shear forcesand externalforceson a unit volume air in the wind form.

259

where r,1 =pu2,1 and r,2 =pu2,9. are the shear stresses under and above the rotor layer, respectively, ct has been defined in general terms previously. Outside the rotor layer it is assumed that ct-0. THe integral in Eq. (22) may be evaluated as h2

Ar= ~ pot u2 dz,~flc't u2 ,

(23)

hi

where C't -~ I CTAr/x2-~ ½CTICr2r/ (Srr)2= CTIC/2S2

(24)

and Uh is the asymptotic, spatial average of wind speed at hub height within the wind farm. CT is the thrust coefficient of the individual wind turbines. By way of the assumption of logarithmic wind profiles over and under the rotor layer the friction velocities of the two profiles are linked by 2 ----ptt.1 2 +pu~c't pU .2

(25)

where u.1 and u.~ are the friction velocities under and over hub height, respectively, and To is the average of the shear stress over and under hub height. For the layer below hub height the roughness length is known, and we get

Uh _ l l n ( h ~ = K,=~u.,_ uh , u.1 k \Zo/ K1

(26)

For the outer layer, the geostrophic drag law is applied:

G / 2 k - u a _ ~1 In ( ~ )

= K 2 = ~ u . 2 ~-

U,2

G/2X--Uh K2

(27)

Inserting Eqs. (26) and (27) in Eq. (25) provides and equation for determination of the spatially averaged wind speed at hub height: K2

J-\-~/

"

This equation may be solved with respect to uh:

G -l+K2(c't+l/K) 2 Uh=-~ 2K2(c~ + 1 / K ~ _ I / K 2 •

(29)

Having determined Uh friction velocities u.1 and u.2 may be calculated. In the real world, it is not likely that the wind profiles are logarithmic all down/up to hub height. We therefore replace the the logarithmic profiles in part of the rotor layer with a linear profil~ in such manner that the resulting profile and its derivative are continuous, and so that the part-profiles joint at heights z = h - Az and z = h + Az, implying that Az = It (u.2- u.1 ) / (u.~ + u.2 ).

260

Figure 4 show winds profiles for different separations of the wind turbines. As described, the profiles are logarithmic at some distance under and over hub height. Investigating the wind speed reduction at hub height, R, = Uh/Uo,where Uo is the free stream wind speed at hub height, it turns out that R,, is little effected by the the magnitude of the geostrophic wind, and that the speed reduction is equally dependent on hub height and roughness parameter. However, realistic hub heights vary at the most a factor 5 while the roughness parameter may vary, say, up to a factor 100. Empirlca!ly, i~ is found that R. may be modelled fairly accurately by the expression

R. ~ InIe_ I

v

(30)

~R+c,t) ,

where 7--O.025/ln(h/zo) 1/3 and c't IECT/s2 and e is the base number of the natural logarithm. Figure 5 compares the empirical expressions, Eq. (30), with the results of the model for terrain roughness %=0.1 m and Zo=0.001 m of the terrain in which the wind turbines are placed. It is seen that the wind speed reduction is considerably larger (R, smaller) for small Zothan for larger values of Zo. =

1 O0

E ¢-

._~ 10 "1-

012345678910 Wind Speed (m/s) Fig. 4. Wind profiles according to the model; Separations: s - 6 , 12, 24, or, hub height h=22 m, wind turbine thrust coefficient Ct-0.88.

261 1.0

/

0.9

//

0.8 ¢-

.o 0.7 Q

o.s

"o 0.5 o

~o.4

ffl

-" Q 0.3 0.2 0.1 0.0

......

ib

.......

'oo

Seporotion s

I • Zo=O.O01 • Zo=0.1

I

Fig. 5. Comparison of model estimated of R,, and the empirical expression Eq. (30). G--15 m, CT = 0.88, h = 22 m. The dotted line is the speed reduction factor emerging when employing Prandl theory for estimating the roughness as it appears at large height~.

The dashed curve is the roughness length obtained employing Prandl theory~ to estimate the roughness from the distributed drag coefficient. Modifying the expression obtained by Prandl theory so that the roughness length approach the terrain surface roughness when the wind turbine spacing increase, the following "outer" roughness length is obtained: h

Z0'2-"exp(~/c~

+-k/ln(zo/h))"

(31)

The expression matchvs exactly the result of the presented model. By introducing the value of Eq. (31) in the equations of the model, an expression more accurate - b u t also more complicated - t h e n Eq. (30) may be obtained.

4.2.Turbulence The structure of turbulence between the wind turbines is complicated and not at all horizontally homogeneous since the nearest individual wind turbine wake upstream of the spatial position considered must be expected to influence the turbulence significantly. 3Assuming the wind profile is logarithmic all down to the terrain surface results in a roughness length as observed from a hrge height of Zo2-exp ( - k / v / ~ t ).

262

Despite this, it may be interesting to probe the spatially averaged properties of the fluctuations in the u-component of turbulence. Following the assumption of logarithmic wind profile under and above the rotor layer, the standard deviation of the along wind component of wind speed (au) is constant outside the rotor layer. As it turns out, the mean value of the standard deviation under the rotor layer aul ~a.u.1, and over the rotor layer, cru2~a.u.1 , is approximately equal to the value of au at hub height outside the wind farm. Therefore, assuming linear variation in c~. over the rotor layer, one could assumed that at hub height a~ is unchanged relative to value outside the wind farm. "['he ratio au2/aul is plotted in Fig. 6 as a function of machine separation. Whereas the average turbulence level at hub height is not effected by the presence of the wind turbines, the difference if turbulence level at the top and bottom of the rotor layer is strongly increasing when wind turbine separation is decreasing. Though somewhat crude, the spatial average scale of turbulence at hub height may be estimated by means of Eq. (21). The model suggests that at hub height (1) the turbulence level is approximately the same as outside the wind farm (au ~ a.,o = a.u., where by approximation the constant a ~ 2.5 and u. is the free flow friction velocity), and (2) the vertical flow shear inside the wind farm and the shear tension approximate the free flow shear: du/dz.~ ( u . / k ) / ( 1 - h), 4.0

1.0

3.5

0.9

0

'O.B '0.7

""D 2.0

'0.6 ~o~ 0.s .~

0

o.4

o 3.0 n, 2.5

7 1.s, 2 ,

,>

10.3

1.0 ~

'0.2

0.5 i 0.0 1 0

0.1 I0

15

~eparation,

I --m-Turbulence R.

20 Rotor Radii

25

0.0 30

--6"-Length Scale R. I

Fig. 6. Ratio of standard deviationof the u-component wind speed under and above the rotor layer (a.Ja.,) and the scale of turbulence at hub height (L} as a function of machine separation. CT----0.8,Cp=0.4, h=22 m and zo=0.01 m.

263

and r~pu 2. Introducing these expressions for f and dividing by a.u., one gets 1/3

3

_

du/dz in Eq. (21) and

xl/3

au,0

L

Akha•,}

\

u~2

~-1

~1,

(33)

where h(ct-Cp) =c~-Cp - ½~(CT-- CD)/s2 and uh/u, =Ruuo/u.=Ruk In(h/ Zo), where CT and CD are the thrust and power coefficients of the individual machines and Uo is the free stream wind speed at hub height. Also, from previous deductions (Eqs. (10) and (11)) we have that the turbulent scale at hub height outside the wind farm is Lo~Akha3.. Introducing these equalities and approximations in Eq. (33), the following expression for the ratio of the scale of turbulence at hub height inside the outside the wind farm is obtaiiDed:

RL_Lo R,,~, 10

s2

IFR,, k2 L I n ( h ) ] 3 + 1}-'=~

CT-CP[R ,_/'h'~']s,,l_, s2 L

"

(34) ,o=, '""'

The normalized scale of turbulence in the wind farm is showr~ in Fig. 6 for realistic values of CT and Cp (in this case not Betz values). It is seen that R~ varies from 1/4 to 1/2 when separation varies from 5-10 rotor radii to 30 rotor radii. Thus, for typical wind turbine separations of 14-16 rotor radii the average scale of turbulence is found to be approximately 1/3 of the scale of turbulence in the free stream, though dependent on the specific values of CT and Cp. Typically, it is found [3,12] that the turbulence scale in the wake of one wind turbine is approximately equal to the rotor diameter, which in turn is approximately the same as hub height or 1/10 of the free stream scale. Data are not at present available to verify the above average considerations consequently it should be stressed that Eq. (35) must be considered as only an indication of the overall turbulence scale at hub height. 5. Comparison with momentum theory In Fig. 7, the asymptotic wind speed reduction factor Ru as determined by the present model is compared with the reduction factor for an infinite, single row of wind turbines, Eq. (5), and the output of a more elaborated computer code based on the same theory, working on a 10 × 10 array of wind turbines. In both cases, the reduction factor has been corrected to account for the fact that wind farm models refer to the wind speed immediately upstream of the wind turbines whereas the model inhere predicts the asymptotic, overall mean wind

264

0.9 _

O

.

B

O.71

I

0"51

I

0"61

I

;°-'1

F/

1

0

0

I o

I

,/ 2

4

NOJ

I

6 8 10 12 14 16 18 20 Separation - Rotor Radii

--m-- BL

-,,at,- PCrnodeB I

Fig. 7. Comparison of model with the single row model, Eq. (5) and the corresponding P C computer code. CT = 0.88, Cp = 0.59, zo = 0.025 m, rr= 1i m and h - 22 m.

speed. (In the shown example with Betz values for CT and Cp, the correction factor becomes I - ~v/2/s.) As seen, the model is in good agreement with the models based on considerations regarding momentum balance in the individual wakes and subsequent superposition of the wakes. One implication of the present model is that the apparent roughness of an infinitelylarge wind farm isnot dependent on the wind direction.This isnot the case for finitesizewind farms, where the array efficiencyactuallymay be stronglydependent on whether wind directionis parallelto the rows of wind turbines or not. In the comparison, Fig.7, ithas been anticipatedthe wind directionbeing parallelto the rows of machines constitutes a situation comparable to the infinitelylarge wind farm. 6. Conclusion

Due to the strong singularities, which the individual machines in a wind farm constitue, it is difficult to generalize the flow description. In the paper, an attempt has been made to evaluate spatially averaged characteristics of the flow for possible use in specification of design criteria. A boundary layer model for the probable asymptotic wind speed in a wind farm has been suggested, in principle only valid in the centre of a very large array of wind turbines. However, experience indicate that output from wind turbines in wind farms stabilizes already 3-4 rows into the array, and conse-

265

quently the model may have practical implications for calibration of models of finite size wind farms. As expression for the apparent roughness of the wind farm has been developed. Considerations concerning turbulence have been presented, indicating that the spatially averaged turbulence level is strongly increasing from the lowest to the highest position of the blade tips. References 1 P.B.S. Lissaman and E.R. Bate, Energy effectiveness of arrays of wind energy conversion systems, Aerovironment Inc. Report AV FR 7058, Pasadena, USA, 1977. 2 I. Katic, J. Hojstmp and N.O. Jensen. A simple model for cluster efficiency; EWEC 86, Rome, Italy, 7-9 October 1986. 3 J.F. Ainslie, Calculating the flowfield in the wake of wind turbines, Advances in Windfarming - Proceedings of the Int. Conf. on Windfarms, Leeuwarden, The Netherlands, October 1316, 1987. 4 N.O. Jensen. A note on wind generator interaction, Rise Report Rise-M-2411, Roskilde, Denmark, November 1983. 5 W. Brutsaert. Comments on surface roughness parameters and the height of dense vegetation, J. Meteoroi. Soc. Jpn., 53(1 ) (1974-75) 96-97. 6 P.S. Jackson. On the displacement height in the logarithmic velocity profile, J. Fluid Mech., 111 (1981) 15-25. 7 F. Engelund. Hydrodynamics - the mechanics of Newtonian fluids (in Danish), Den Private Ingenicrfond, the Technical University of Denmark, 1968. 8 J.C. Kaimal, J.C. Wyngaard, Y. Izumi and O.R. Cote. Spectral characteristics of surface layer turbulence, Quart. J. Roy. Meteorol. Soc., 98 (1972) 653-689. 9 R.M. Cionco. A wind-profile index for canopy flow, Boundary Layer Metereol., 3 (1972) 255263. 10 A.S. Thorn. Momentum absorption by vegetation, Quart. J. Roy. Soc., 97 (1971) 412-428. 11 R.M. Cionco. A mathematical model for air flow in a vegetative canopy, J. Appl. Meteorol., (Aug. 1965). 12 J. Hsjstrup. Measurements of structural loads in the Ncrrekmr Enge wind farm, Proc. of the BWEA CGi,ference, Swansea, UK, April 1991.