Wind turbine with a wind velocity measurement system

[0002] This application is the national phase under 35 U.S.C. §371 of PCT International Application Ser. No. PCT/DK98/00125 which has an International filing date of Mar. 25, 1998, which designated the United States of America.

[0003] The present invention relates to a wind turbine comprising a wind velocity measurement system for determination of air velocities in front of the wind turbine.

[0004] Wind turbines are becoming an important source of electric energy. The majority of wind turbines are connected to the electric grid. As the number of larger turbines increases the need for reliable local wind measurements for power curve determination and for advanced control becomes more important.

[0005] Having a history of millennia and a century of obscurity following the introduction of fossil fuel engines, the wind power technology had its revival in the mid-seventies as a consequence of the oil crises. Various traditional applications such as water pumping and grain-grinding and small wind turbines for charging of batteries in remote areas have, in modernised form had their renaissance. However, a large fraction of the global energy supply is now delivered in the form of electricity and thus the major effort in wind power technology was placed in the development of electricity generating machines. Therefore, nearly the entire manufactured wind turbine capacity is now so-called grid connected wind turbines.

[0006] The modern “stock-average” wind turbine has an aerodynamic rotor with two or—most frequently—three blades. The rotational speed of the rotor is one-half to one revolution per second; through a gearbox the speed is increased to a level suitable to an electrical generator which is typically of the induction type: the rotor of the (electric) generator is magnetised by the grid and once the connection to the grid has been established the rotational speed is kept constant within 1-2% over the range of operational wind speeds, 3-25 m/s. The nacelle (machine cabin) with the gearbox and electric generator is mounted on a tubular or lattice tower through a motor-driven yaw system which by means of a wind direction sensor keeps the wind turbine facing the wind. When the blades are fixed to the rotor hub the machine is said to be stall-regulated because the blade profile angle to the wind has been selected so that the blades stall and loose power when the rated power of the electric generator is reached. If the blades can be pitched the machine is said to be pitch-regulated: when rated power of the electric generator is reached, the blades are pitched to limit the aerodynamic power.

[0007] Wind turbines operating with variable speed are also known.

[0008] The grid connected wind turbine is equipped with a control system providing smooth cut-in and cut-out of the grid, correct yaw positioning of the rotor, and pitch setting of the blades, stopping of the machine if voltage is too high or too low or if loads on the phases are not balanced, stopping of the machine if unacceptable vibrations are experienced or the temperature of gearbox or generator is too high etc. A number of sensors are connected to the control system; sensors monitoring the air flow may be a wind vane and a cup anemometer mounted on top of the nacelle—just behind the rotor and the flow is therefore rather disturbed compared to free-flow conditions.

[0009] The grid-connected wind turbines built and sold commercially in the seventies had installed capacities of 20 to 50 kW and aerodynamic rotor diameters of 5 to 15 m. In 1997, the typical commercially available machines have installed capacities of 500 to 800 kW and rotor diameters of 40 to 50 m. The next generation of machines presently being tested are 1 to 1.5 MW. This urge to increase size probably stems from belief in “economy of scale”, i.e. the larger the cheaper, and from the fact that it is easier to fit into the landscape one 1 MW machine than ten 100 kW machines.

[0010] Cost of energy from wind power plants depends foremost on the wind energy resources, i.e. the kinetic energy flowing through the area covered by the wind turbine's rotor (per sec.): $\begin{array}{cc}{P}_{\mathrm{total}}=A_r{\int }_0^{\infty }1/2\rho u^{3}f\left(u\right)du,& \left(1\right)\end{array}$

[0011] where A_ris the area swept by the rotor, ρ is the air density, u is the wind speed and f is the frequency distribution of wind speed. The consequence of the available energy being a function of the cubed wind speed is that if the annual mean wind speed at one site is 20% higher than at another site the available energy is {(1.2)³−1}≅70% higher. Thus, moving from a site with an annual mean wind speed of 5 m/s to a site with 6 m/s the cost of energy (potentially) drops 1/1.7≅58% thus, illustrating the vast importance of the wind climate to the economy of a wind power plant.

[0012] The other factors are cost of the wind turbines, operation and maintenance and the efficiency of the wind turbines. Over 20 years, improved aerodynamic and structural design tools and production methods have reduced cost of energy—disregarding the wind climate—by 70%. By 1997 the cost of energy at European sites with high wind speeds was 0.03-0.04 USD/kWh, and at “medium-good” wind sites 0.05-0.07 USD/kWh. Reduction in cost of energy becomes increasingly difficult. However, it is expected that the next 30 years will bring a further reduction of 20-30%.

[0013] FIG. 1 shows a conceptual lay-out for a wind turbine with control based on upstream wind measurements,

[0014] FIG. 2 shows a normalised output power, Cp, (lower traces) for a wind urbine operating at fixed and variable rotational speed, respectively (upper traces),

[0015] FIG. 3 shows the additional energy due to turbulence,

[0016] FIG. 4 shows a basic closed-loop control system,

[0017] FIG. 5 shows a contour plot of shaft torqueas a function of wind speed and pitch angle (WD34°35 rpm),

[0018] FIG. 6 shows an overview of a wind turbine from a control point of view,

[0019] FIG. 7 shows a nomalized output curve, Cp, curve for WD34 for a pitch angle of 1 deg,

[0020] FIG. 8 shows an open loop estimation of wind speed,

[0021] FIG. 9 shows a closed loop estimation of wind speed,

[0022] FIG. 10 is a blocked schematic of a laser anemometer system,

[0023] FIG. 11 shows the number of photons/sec scattered by particles in a single spatial mode versus the optical wavelength assuming that 1 watt of light power is transmitted in a counter propagating mode otherwise identical to the light collection mode. The parameter v is a characteristic parameter of the Junge particle size distribution. A large v implies a relatively large amount of small particles. Note that the number of photons per second peaks at around 10 μm almost independently of v. The wavelength dependence of the refractive index has not been incorporated,

[0024] FIG. 12 shows a layout of the optics of a back-scattering reference beam configuration. The polarising beam-splitter in conjunction with the wave plates ensures efficient utilisation of the available laser power,

[0025] FIG. 13 shows a layout of an autodyne or self-mixing configuration. The scattered light propagating in the same mode as the transmitted light—except for the direction of propagation—enters the laser and perturbs the internal power level. This is detected by the photodetector at the rear end of the laser. It is noted that a it is very important to avoid reflections that would enter the laser in a mode that could be amplified in the laser,

[0026] FIG. 14 shows an autodyne system with a tapered waveguide laser,

[0027] FIG. 15 shows a beam expander based on reflective optics. The location of the measuring volume in a plane perpendicular to the initial optical axis can be changed by tilting the largest of the two mirrors. The distance to the measuring volume can be changed by a displacement of the large mirror. It is important that it is the large mirror that is moved.

[0028] FIG. 16 shows a signal processing layout. The frequency locked loop (inside the dashed line) may be implemented in the form of an application specific integrated circuit. The microprocessor sets the range as well as the loop speed and will also control the initiation of the lock. The input filter is envisaged to be based on traditional analogue circuitry. It serves to eliminate the dc signal from the detector and to confine the noise bandwidth,

[0029] FIG. 17 shows the beam parameters: d_tis the diameter of the transmitter/receiver aperture, d_ois the focal beam diameter, L is the distance to from the anemometer to the centre of the measuring volume and 2Δz is the length of the measuring volume,

[0030] FIG. 18 shows a simulated detector signal assuming 10% turbulence, The third line concerns the development of advanced control strategies and systems that will improve efficiency, reduce static as well as dynamic structural loads, reduce noise, and facilitate the grid connection in terms of power quality. Also the application of power electronics, allowing for variable-speed operation, will probably be common in the future,

[0031] FIG. 19 shows characteristic parameters in conjunction with velocity measurement from a wind turbine,

[0032] FIG. 20 shows plots of H_R(k) for various weighting functions across the rotor,

[0033] FIG. 21 shows a plot of the correlation function σ(Xi₁),

[0034] FIG. 22 shows a plot of the correlation function σ⁰(Xαi₂), and

[0035] FIG. 23 shows a plot of the error function erf(x).

[0036] It is an object of the invention to provide a wind turbine with improved efficiency, reduced static as well as dynamic structural loads, reduced noise, and facilitating grid connection with improved power quality.

[0037] According to the present invention the above and other objects are fulfilled by a wind turbine comprising a wind velocity measurement system for determination of air velocities in front of the wind turbine and having an output for provision of determined air velocity values.

[0038] The wind velocity measurement system may comprise a laser anemometry system having a laser for emission of a laser light beam for illumination of particles in the air, a detector for detection of light emitted from particles that have interacted with the illuminating laser light, and a signal processor connected to the detector for generating a velocity signal corresponding to the velocity of the particles based on the detected light.

[0039] According to a preferred embodiment of the invention, a wind turbine is provided comprising adjustable pitch blades, a pitch controller for setting blade pitch and having an input, the blade pitch being set according to a signal provided at the input, and a control system for controlling the output power generated by the wind turbine by control of blade pitch. The control system may have a control output that is connected to the input of the pitch controller, a velocity input that is connected to the output of the velocity measurement system, and a pitch processor for calculation of the desired blade pitch based on determined air velocity values and generation of a corresponding pitch control signal that is provided at the control output.

[0040] It is an important advantage of the invention that the determination of air velocities in front of the wind turbine compensating for delay in blade pitch setting so that the actual blade pitch setting corresponds to the actual wind velocity at the blades during operation of the wind turbine.

[0041] According to another embodiment of the invention, a method of wind turbine power curve determination is provided comprising the steps of mounting on the wind turbine a wind velocity measurement system for determination of air velocities in front of the wind turbine, connecting a power measurement system for determination of the output power generated by the wind turbine between a power generator of the wind turbine and a load receiving the power generated by the wind turbine, and recording corresponding values of wind velocity from the wind velocity measurement system and output power from the power measurement system.

[0042] According to still another embodiment of the invention, a method of controlling output power from a wind turbine is provided comprising the steps of mounting on the wind turbine a wind velocity measurement system for determination of air velocities in front of the wind turbine and providing determined velocity values to a power control system of the wind turbine.

[0043] Optimisation of wind turbine efficiency—and in turn measurement of the efficiency—includes accurate measurement of wind speed in the immediate vicinity of the wind turbine. A laser anemometer mounted on the wind turbine could serve both purposes possibly with better accuracy than the presently used cup anemometer, see FIG. 1.

[0044] Wind turbines may have pitchable blades and operating with variable rotational speed of the rotor providing optimal energy extraction and improvement of power quality. Under given flow conditions—mean wind speed, turbulence shear etc. —there is an optimal setting of rotor speed and blade pitch. However, the response time of the wind turbine is several seconds and thus, a short-term prognoses of the flow characteristics is needed. The laser anemometer may potentially deliver such forecast. The characteristics of a wind turbine operating at a constant and variable rotational speed, respectively, is illustrated in FIG. 2. The term C_pis a normalised power output of the wind turbine, defined by the following expression: $\begin{array}{cc}{C}_p\left(u\right)=\frac{p\left(u\right)}{\frac{1}{2}\rho u^3{A}_r}& \left(2\right)\end{array}$

[0045] where p(u) is the so-called power curve of the wind turbine: the power out as a function of wind speed. The efficiency of fixed-speed wind turbine is a function of wind speed and is maximum only in a narrow wind speed band.

[0046] Variable-speed operation will expand the range of high (maximum) efficiency, and the laser anemometer facilitates optimal control of speed with regard to efficiency—and potentially also for structural load and noise reduction.

[0047] Further, the laser anemometer facilitates (retrospective) assessment of production capability of wind power plants, being important to the credibility of the wind power technology in the financial world. The measurement of wind turbine efficiency involves measurement of wind speed at approximately 50 m above the ground and 100-200 m upstream, where the flow is undisturbed of the machine itself—FIG. 1. The laser anemometer may measure the wind speed without having to erect a separate tower; in addition, the laser anemometer may provide a spatial average of wind speed and turbulence accurately representing forces on the wind turbine. This would contribute to an improvement of the accuracy of performance assessment of wind power plants.

[0048] Below, characteristics are outlined of the boundary layer flow field. Such characteristics are assumed important in determining the performance of a wind turbine and in selecting a strategy (for maximum energy extraction) for the rotational speed of rotor and/or blade pitch angel.

[0049] Power output from a wind turbine may be predicted based on knowledge of the wind speed frequency distribution at the wind turbine site and of the so-called power curve of the wind turbine. The frequency distribution shall be given for the wind turbine's hub height, and the power curve is measured with an anemometer elevated to hub height upstream of the wind turbine. In both cases 10 min averages of wind speed (and power) are applied. As the same averaging time is used for both power curve and frequency distribution, the end-result is independent of averaging time. However, in the real world turbulence, flow shear, air density, etc, may be different for the site where the power curve was measured and the site on which the wind turbine is to be positioned. Experience indicates that in fact these parameters are important and uncertainties are introduced when disregarding them.

[0050] As stated, when measuring power curves turbulence is disregarded as being of lesser importance. However, in stationary flow with mean U and (for illustration) isotropic turbulence (the flow has identical characteristics in all directions, thus also the same fluctuations in the three directions: σu=σv=σw) the kinetic energy in 1 m³of air is $\begin{array}{cc}E=\frac{1}{2}\stackrel{\_}{\rho {\left(U+u\right)}^2+v^2+w^2}=\frac{1}{2}\rho U^2+\frac{3}{2}\rho {\sigma }_u^2& \left(3\right)\end{array}$

[0051] where overbar denotes average (note that it is assumed that the air density is constant); the energy flux (power) through an area of 1 m²perpendicular to the mean wind direction is

[0052] P={overscore ((U+u)E)}=½ρ{overscore ([(U+u)²+L =v²+L +w²+L ](U+u))}≈½ρU³(1+5I_u²),  (4)

[0053] where I_u=σu/U is the turbulence intensity. Thus, for turbulence intensities 10-20% the contribution to “per-m²-energy-flux” from turbulence is 5-20%. More realistically, in the atmosphere, the coefficient to I_uis of the order 3-4 rather than 5 since σu is larger than the other two components. Thus, for turbulence intensities 10-20% turbulence increases the energy flux per m²with 5-15% additional to what is calculated from the mean wind speed alone. Disregarding lateral and vertical components of turbulence the energy flux becomes

[0054] P=½ρ{overscore ((U+u)²+L (U+u))}≈½ρU³(1+3I_u²)  (5)

[0055] Here, the contribution to energy flux from turbulence of 10-20% is 3-12%. Performance data indicate that power output may vary up to 8% for variation in turbulence intensity from under 5% to over 15%. In performance standards, it is presently not taken into account that the energy flux depends on turbulence intensity.

[0056] Thus, improvement of the accuracy of performance testing may amongst other things be improved by also measuring the turbulence that the machine is exposed to. To do that a suitably spatial/time average of turbulence must be defined and measured. The laser anemometer facilitates contribute measurement of turbulence.

[0057] Evidence points to higher power output when turbulence increases even for fixed-speed machines. If the rotor speed can be varied a larger fraction of the turbulent energy may be extracted, and it may be assumed that a continuous prediction of the instant wind speed—such as supplied by a laser anemometer—may help select the optimal rotor speed at any given time.

[0058] Considering the potential improvement of energy extraction by variable-speed operation of the wind turbine, it is useful to investigate the potential gain. The gain may be divided into what is gained by a relative slow regulation (time constant of the order of 10 min) and a faster regulation (time constant of a few seconds).

[0059] Usually, climatology wind statistics is based on 10 minute average values of the wind speed. It is interesting how the averaging time influence the estimation of energy flux. Considering only the U-component of the wind field, the annual average per-m2 energy flux is $\begin{array}{cc}P=F_E\left(k\right)A^3=1/2\rho \Gamma \left(1+\frac{3}{k}\right)8.76A^3=1/2\rho 8.76\frac{\Gamma \left(1+3/k\right)}{{\Gamma }^3\left(1+1/k\right)}{\stackrel{\_}{U}}^3\left[\mathrm{kWh}{m}^{-2}{y}^{-1}\right]& \left(6\right)\end{array}$

[0060] where A and k are the Weibull parameters and F_Ea function tabulated in Risø Report 1981 by Petersen et al. If the turbulence intensity (I_u=σ_u/U) is zero, the energy flux is given by eq. 6 applying values of A and k determined from measured 10 minute averages of wind speed. If the turbulence intensity is not zero, the annual energy flux will be under-estimated. The difference in Weibull parameters between 10 minute averages and instantaneous values of wind speed has been evaluated. The relative differences, δA/A and δk/k, respectively, are for small values $\begin{array}{cc}\frac{\delta k}{k}=-\frac{k{I}_u^2}{2\left(\psi \left(1+2/k\right)-\psi \left(1+1/k\right)\right)},\bigwedge \frac{\delta A}{A}=\frac{1}{k}\psi \left(1+\frac{1}{k}\right)\frac{\delta k}{k}& \left(7\right)\end{array}$

[0061] where Ψ is the digamma function: $\begin{array}{cc}\psi \left(x\right)=\frac{d\ln \Gamma \left(x\right)}{dx}& \left(8\right)\end{array}$

[0062] It is seen that δA and δk are functions of k and Iu. The energy flux neglected when using 10 minute average values is $\begin{array}{cc}\frac{\delta P}{P}=\frac{F_E\left(k+\delta k\right)·A^3\left(k+\delta C\right)-F_E\left(k\right)·A^3\left(k\right)}{F_E\left(k\right)·A^3\left(k\right)}& \left(9\right)\end{array}$

[0063] For the value k=1.5 and 3, the eq. 9 is plotted in FIG. 3 together with the expression of eq. 5.

[0064] The points are the increase calculated by means of eq. 9 and the solid line is the simple eq. 5. It is seen that the energy flux increase is only slightly dependent on k and that the eq. 5 gives a good estimate not only in stationary flow but also in the atmospheric boundary layer, of the fraction of energy flux left out when using 10 minute average as a climatological base number.

[0065] At realistic turbulence levels of 13 to 18% the energy flux ignored by using 10 minute statistics to calculate frequency distribution parameters is 5 to 10%.

[0066] Since the turbulent kinetic energy usually is ignored, the non-constant turbulent kinetic energy appears as an intrinsic uncertainty, e.g. an uncertainty which is not eliminated by improvement of test equipment or duration of test.

[0067] Thus, in highly turbulent air flows, the energy flux can not be accurately determined from the mean wind speed only. However, there are other intrinsic contributions to be considered: the wind turbine extract energy from an extended area (the swept area of the rotor), the single-point measurement at hub height may not sufficiently well represent the flux of energy experienced by the rotor. Thus, also horizontal and vertical wind shear and misalignment (yaw error) may cause errors in the power curve measurement. It has been found that shear and yaw error may cause variations in predicted annual energy output of more than 5%.

[0068] The largest, present-day wind turbines have two or three blades of lengths up to 40-45 m; the rotational speed of the blades is high, with blade tip speeds of the order 60 m/s. The blade profiles are similar to airplane wing profiles, the purpose of which is to reduce drag and in turn energy conversion losses.

[0069] For the ideal rotor—e.i. a rotor with no drag losses and with infinite high rotational speed—operating in a non-turbulent flow, it can be shown that the maximum of the rotor efficiency defined as $\begin{array}{cc}{C}_p=\frac{P}{1/2\rho u^3{A}_r}& \left(10\right)\end{array}$

[0070] (where P is extracted power, A_rthe swept rotor area, ρ air density and u wind speed) is $\begin{array}{cc}\max \left[C_p\right]=\frac{16}{29}\approx 0.59& \left(11\right)\end{array}$

[0071] Further, it can be shown (still employing simplest possible theory) that the maximum efficiency is obtained when the flow speed behind the rotor is reduced to ⅓ of the upstream flow speed and to ⅔ at the rotor disc. Betz theory is inadequate when dealing with real rotors; however, it does point to the important fact that not only is the flow as measured by the presently used nacelle cup anemometer disturbed by the nacelle itself but the flow speed is inherently linked to the performance of the rotor: e.g the higher performance the lesser wind speed is measured on the nacelle top, implying that the value of this instrument for performance measurements is limited.

[0072] Presently, the so-called blade element theory is the method applied for aerodynamic load -and power calculations, though with various auxiliary methods to take into account e.g. turbulence, stall behaviour of the blades and the induction or wake delay. The blade element theory quite accurately models power output and loads under stationary and quasi-stationary conditions when the flow around the airfoils are far from stall. The blade element theory does not account for dynamic flow behaviour. Therefore, to model the time-wise development in aerodynamic power, modifications have to be applied.

[0073] Below, the conventional method for power performance measurements is reviewed and intrinsic uncertainties which may be reduced or eliminated by employment of nacelle mounted laser anemometry are pointed out.

[0074] Since the IEA recommendation for power performance testing was issued in the late 1970' ies it has been found that actually unacceptable large uncertainties were attached to the tests: even under fairly ideal test conditions is was estimated that the uncertainty could reach a level of 10 to 20%.

[0075] Due to lack of a genuine commercial market this fact was not seriously acknowledged until a rapidly increasing commercialisation has hit the wind power industry in the resent few years. In reality, firm contracts between manufacturers and developers now implicitly call for uncertainty levels for performance testing as low as 1-2% in terms of standard uncertainty, which is unrealistic—even though the understanding of the problems has improved.

[0076] As of now a “Committee Draft” for power performance testing, IEC 96, has been prepared. The document has not yet obtained the status as an international standard but will probably do so in the near future. The document offers a procedure to be followed when measuring power curves, as well as a method for assessing uncertainties. The method and the proposed uncertainty analysis are outlined below.

[0077] While the measured power curve may be used for other purposes, the IEC document is directed toward application to calculation of annual energy production (AEP). The power curve is basically established by collecting pairs of wind speed and power measurements, though also additional measurements are made.

[0078] Test site: shall be assessed for sources of flow distortion, mean terrain or other features causing differences in wind speed at the position of the wind turbine and the met tower. At least topographical variations, other wind turbines, and obstacles shall be taken into account in the analysis.

[0079] Meteorological tower: shall be placed between 2 and 4 rotor diameters (D) from the wind turbine; more than 2D to avoid stemming disturbances, and less than 4D to maintain good correlation. 2.5D is recommended. Measurements are acceptable for approx. 130° around the case where the wind is from met tower to wind turbine.

[0080] Correction factor due to flow distortion: if the test terrain adheres to specified requirements no corrections need to be made and the uncertainty related to flow distortion shall be set to not less than 2-3% depending on separation between met tower and wind turbine. If the test terrain does not adhere to the specified requirements of an ideal test terrain site calibration shall be performed. If a computer model is applied the uncertainty shall be assumed to be ½ of the correction made to wind speed; if measurements used for site calibration the uncertainty is assumed to be ⅓ of the correction.

[0081] Electric power: The sensor(s) to measure net electric power shall be of class 0.5 or better.

[0082] Wind speed: Shall be measured with a cup anemometer mounted at the wind turbine's hub height H_h“2.5% above the ground in the met tower. It is specified how to mount the instrument so as to avoid flow distortion from the met tower itself. The distance constant of the anemometer shall be less than 5 m/s, the instrument shall be calibrated and the uncertainty stated.

[0083] Wind direction: Shall be mounted at H_h“10% in the met tower and the accuracy shall be better than 5%. The purpose of the instrument is to single out data, where the met tower is in the wake of the wind turbine.

[0084] Air density: Should be (indirectly) measured at H_hat the met tower, by measurement of pressure and temperature in accordance with ISO 2533.

[0085] Precipitation: Shall be monitored to distinguish between wet and dry periods.

[0086] Wind turbine status: A parameter indicating whether the wind turbine is ready for operation or not.

[0087] Data acquisition system: Sample frequency shall be at least 0.5 Hz. It shall be ensured that uncertainties related to the data acquisition system are negligible compared to other uncertainties.

[0088] Measurement procedure: During the measurement the wind turbine shall be in operation and the above parameters shall be sampled at 0.5 Hz or faster. From the data mean values of 10 min periods of consecutive data sets, standard deviation, maximum and minimum shall be computed and stored. At least 180 hours of data shall be collected, and there must be at least 30 min of data in each ½ m/s wind speed bin.

[0089] Data correction: The data of the data base must be corrected (calibration errors, site calibration correction etc.)

[0090] Data normalisation: Data shall be normalised to two reference air densities: one rounded off to the nearest 0.05 kg/m³and one to the ISO standard atmosphere 1.225 kg/m³; the precise method is given in the IEC document.

[0091] Power curve: The corrected and normalised data are binned in 2 m/s wind speed bins and averages of wind speed and power are formed in each bin. These sets of bin-averaged values of wind speed and power constitute the power curve.

[0092] Annual energy production (AEP): The AEP is estimated by applying the measured power curve to annual mean wind speed, assuming the 10 min average wind speeds to be Rayleigh distributed. The AEP is what in metrology is termed the measurand.

[0093] Knowledge—or at least awareness—of the component uncertainties and their mutual relations in an experiment is the key to improvement of the test.

[0094] A guide to uncertainty analysis, ISO 1993, was issued recently, having all major standard organisations behind it. One of the main points of the guide is to disregard the concept of bias uncertainty and replace it with “correlated uncertainty”. Doing so, all uncertainties may be treated equally, e.g. added in the same way, which is a major improvement. However, it somewhat complicates the presentation of the issue.

[0095] The true value (the concept of “true value” is also estranged in ISO 1993; it is used in the present disclosure anyway) of the measurand is never exactly known, and therefore determination of the uncertainty is a probabilistic problem, i.e. what is the probability that the true value is found inside the interval [μ−σ;μ+σ], where μ is the measured value and σ is the uncertainty measure? In some cases objective measures of the sought probability exist, i.e. when a stochastic variable is measured (with no instrumental uncertainty) and the measured values are independent and form a Normal Distribution with a mean μ and a standard deviation σ, the probability of each measured value being in the interval [μ−σ;μ+σ] is approximately 68%, the mean and standard deviation being estimated from the measurements. In other cases the shape of the distribution and/or the standard deviation may not be explicitly known.

[0096] If the distribution is not Normal, the probability attached to the standard deviation may be different. However, in this annex, it is assumed that all uncertainty components are sufficiently accurately characterised by a standard deviation. Further, it is trivial—though important—to note that the standard uncertainty for e.g. a random variable is not identical to standard deviation of the process, σ. The uncertainty standard deviation is the estimated standard deviation of the final average (binwise or else-wise) of a component in contrast to the standard deviation, σ, of the parameter.

[0097] There are two types of uncertainties. It is assumed that they are only distinguished by the way their magnitudes are determined.

[0098] The first type of uncertainty, category A, is characterised by firm knowledge of the uncertainty distribution, e.g. through the data base of the power performance test or tests made prior to the performance test.

[0099] The distribution of the measured quantity is assumed Normal and the standard uncertainty is estimated by the expression $\begin{array}{cc}u=\frac{\sigma }{\sqrt{N}}& \left(12\right)\end{array}$

[0100] where σ is the standard deviation of the measured process and N is the number of observations. Uncertainties of category A are in principle independent of the person who makes the uncertainty estimate, and the uncertainty of the estimate of the mean value can in principle be reduced to any small quantity by repeated measurements.

[0101] The second type of uncertainty, Category B, is characterized by no or imprecise knowledge of the uncertainty distribution. Category B uncertainties is estimated by other means than measurements. The uncertainty type would typically be associated with “bias” uncertainty; however, if a bias/one-sided error is identified the data should be corrected. What is left after correction may or may not be biased; in any case the sign is not known.

[0102] The (important part of the) uncertainties of instruments is often of Category B: little or no specific information may be available—or a statement of uncertainty limits, for instance “U, is given in data sheets for the instrument. For instruments that are only specified by uncertainty limits, an uncertainty distribution as well as uncertainty standard deviation shall be estimated. Often, a rectangular distribution is chosen. Uncertainties of category B are often based on personal judgment, and consequently they are often dependent on the person who makes the uncertainty analysis. Uncertainties of category B can not be minimized by repeated measurements.

[0103] The square of the standard uncertainty of annual energy production (AEP) can in its general form be expressed as $\begin{array}{cc}{u}_{\mathrm{AEP}}^2=N_h^2\sum _{i=1}^N\sum _{j=1}^N\sum _{k=1}^M\sum _{l=1}^M{f}_i{c}_{k,i}{u}_{k,i}{f}_j{c}_{l,j}{u}_{l,j}{\rho }_{k,l,i,j}& \left(13\right)\end{array}$

[0104] where

[0105] f_i=relative occurence of wind speed in wind speed bin i,

[0106] u_k,i=standard uncertainty of component k in bin i,

[0107] c_k,i=sensitivity factor of component k in bin i,

[0108] M=number of uncertainty components,

[0109] N=number of bins,

[0110] N_h=8760, number of hours per year,

[0111] ρ_k,l,i,j=correlation coefficient between uncertainty component k in bin i and uncertainty component l in bin j.

[0112] For all uncertainty components it is assumed in the IEC document that

[0113] different uncertainty components are uncorrelated, ρ_k,l,i,j=0 for k≠l.

[0114] In addition, for type A uncertainties it is assumed that

[0115] Uncertainty components in different bins are uncorrelated, ρ_k,l,i,j=0 for i≠j.

[0116] This implies that the square of the integrated type A uncertainty can be written as $\begin{array}{cc}{u}_A^2=N_h\sum _{i=1}^N\sum _{k=1}^{M_A}{f}_i^2{c}_{kA}^2{u}_{i,k}^2& \left(14\right)\end{array}$

[0117] where M_Ais the number of uncertainty components of type A and _Au²_i,kis the uncertainty of type A, component k in bin no. i.

[0118] For category B uncertainties, the additional—possibly conservative—assumption is that bin-wise uncertainties of the same component is fully correlated throughout all bins, ρ_k,k,i,j=1.

[0119] This implies that the square of the integrated type B uncertainty can be written as $\begin{array}{cc}{u}_B^2=N_h\sum _{k=1}^{M_B}{\left\{\sum _{i=1}^N{f}_i{c}_{kB}{u}_{i,k}\right\}}^2& \left(15\right)\end{array}$

[0120] where M_Bis the number of uncertainty components of type B and _Bu²_i,kis the uncertainty of type A, component k in bin no. i. The motivation for the assumption of the type B uncertainties each being fully correlated over all bins is that a typical type B uncertainty is the uncertainty of an instrument's calibration constant, which presumably remains constant during the test. The total standard uncertainty on annual energy production is $\begin{array}{cc}{u}_{\mathrm{AEP}}=\sqrt{u_A^2+u_B^2}& \left(16\right)\end{array}$

[0121] In performance testing as in many other tests the type A uncertainties are less important since these are reduced any small number if the duration of the test is long enough, Eq. 12. However, the type B uncertainties remain constant once test method, equipment, site etc. have been chosen. Therefore, improvement of the accuracy of performance testing primarily implies reduction of the type B uncertainties.

[0122] Below, the uncertainties constituting important error sources are listed and improvements of the performance test provided by the laser anemometer are assessed. What presently appear to be the contributing uncertainties in a performance test is outlined in the table below. The list of uncertainty components may well be incomplete, however it reflects the present knowledge. The characteristics of main groups of uncertainties are outlined in the notes to the table. Typical ranges of uncertainties for each uncertainty component are listed in the table.

[0123] The later anemometer appears to be capable of providing improvement options to all uncertainty components, except “obstacles”:

[0124] Instrumental (anemometer): There are indications that calibration accuracy of (cup) anemometers in wind tunnels at present cannot be made better than 1-2% due to blockage effect and lack of natural turbulence. The cup anemometer measures the horizontal speed, S={square root over (U²+L +V²+L )}. While in the natural turbulence both U and V are non-zero, in the wind tunnel V=0. Therefore, the conditions of the wind tunnel are intrinsically different from the real world, which introduces an error. The laser anemometer is a “first-principle” instrument, which in principle does not need calibration, and it will basically measure only the U-component of wind speed, presumably the component experienced by the cosine-responding wind turbine propeller.

[0125] A wind turbine mounted laser anemometer may facilitate experimental site calibration, both in terms of accuracy and in amount of data gathered within a certain period of time. To apply the laser technology for site calibration the instrument must be capable of measuring (at least) the axial wind speed at selected distances from the wind turbine. Also, with the availability of an operational laser anemometer, site calibration will not be a separate test but an integrated part of the performance test.

[0126] Assessing the production capability of a wind farm the power curves are measured on a sample of wind turbines, typically one or two. With the achievable accuracy of the individual power curves, larger samples are needed. However, the per-unit present costs of power curves are high (due to met towers and instruments). A low-cost laser anemometer would make it possible to performance test all units, and thus removing sample uncertainty completely.

[0127] Amount of data: The laser anemometer may collect data from a larger range of wind directions and therefore reduce the integrated measuring time or decrease the random (category A) uncertainty.

[0128] The benefits of an operational laser anemometer will be

[0129] Reduction in direct costs of performance tests.

[0130] Reduction of the uncertainty of the individual power curve as well as of the entire wind farm.

[0131] Lesser financial costs following the lesser uncertainty of the power curve (investment in wind power plants becomes less risky when output is better predicted; therefore a lower interest rate is asked by investors). With a laser anemometer as standard equipment on (large) wind turbines the logistics of performance (and other) tests will be minimised

[0132] (1) Uncertainties directly related to calibration of instruments and their use; Most measurements consist of several measurements in combination.

[0133] (2) Refers to the applied method itself (hub height measurements of wind speed in front of rotor etc.); the potential pitfalls are to some extend dealt with in section 3.1.

[0134] (3) The present method, IEC 1996, highly relies on methods to re-calculate the wind speed measured at met-tower to wind speed at the wind turbine; with the IEC method site calibration uncertainties are potentially severe though very difficult to quantify.

[0135] (4) Analyses, Frandsen 1994, indicate that obstacles (of sizes of the wind turbine) as far away as 1 km may influence the accuracy of the performance test.

[0136] (5) When evaluating plants of a larger number of machines, present practice is to test one or two turbines, which leaves sample uncertainty.

[0137] (6) For both power variability (power curve measurement) and variability when performing site calibration there are (measurable) uncertainties Table 3.1 Summary of uncertainties in power performance testing. Only uncertainties related to “amount of” data are quantifiable. The indicated uncertainty levels are indicative.

[0138] The operation of a wind turbine can be divided in four modes depending on the wind speed:

[0139] Below cut-in

[0140] cut-in to rated

[0141] rated to cut-out

[0142] above cut-out.

[0143] The behaviour of the wind turbine in the two power production modes depends on its type, whether it is stall-, pitch or pitch/variable speed controlled. Considered here are only pitch and combined pitch and variable speed controlled machines. Both these two kinds of wind turbines apply active control mechanisms for limiting the output power. The same mechanisms can be utilised for production optimisation i.e. operating the wind turbine at maximum Cp.

[0144] The basic objective of a control system is to reach a given value of an output signal through manipulation of an input signal. There are two separate ways of controlling a given system: open and closed loop. In the open loop approach the input variable is assigned a value based on a priori knowledge of static characteristics of the controlled system e.g. to have 500 kW from a given wind turbine the wind speed is measured to 15 m/s and then the pitch angle should be 5 deg. and the rotational speed should be 30 rpm based on pre-calculated Cp curves. The actual output is not measured. In the closed loop approach the pitch angle and rotor speed are determined on the basis of measurement of the output power. This ensures that deviations from the assumed characteristics and disturbances are suppressed.

[0145] The closed loop performance is limited by several factors. It is seen, in FIG. 4, that if the phase shift through the loop is 180° and the loop gain is above 1 the fed back signal is actually added to the reference signal and the resulting system is therefore unstable. Also in the figure it is shown how a disturbance, d, can be fed forward, through H₂, in order compensate for the disturbance before it has changed the output through G₂and G₁. It is of course necessary to be able to measure the disturbance if it is to be fed forward.

[0146] The output power from a wind turbine has to be limited in order to limit the loads on the various components. The output can be limited e.g. by applying variable pitch and/or variable speed.

[0147] For pitch control the basic control scheme can illustrated as in FIG. 5. The objective of the control system is to keep the torque constant (it is a, constant speed machine), in the current case to keep the torque at 110 kNm for all wind speeds where it is possible.

[0148] The setting of the pitch angle is controlled by a controller, typically a PI controller with electrical power as input. The electrical power is used for several reasons. It is among other things the power that has to be limited and it is simple, robust and cheap to measure. From a control point of view the electrical power is a more problematic signal to use. The basic problem is that instant torque on the rotor is determined by the instant pitch angle and the instant wind speed as experienced by the rotor and the power the is used for the control loop is the electrical power that appears after the torque has been converted by the gearbox and the generator. The controller then has to calculate the pitch angle that compensates for the deviation in power production and the pitch system has to move to the desired pitch angle in order to limit the power. Unfortunately the very high sensitivity to wind speed changes and the phase lag resulting from the conversion of wind power electrical power makes it very difficult to obtain a power limiting performance that is adequate.

[0149] One of the classical and best ways of improving a control system is to feed forward the signal that causes the disturbance in order to compensate for the effect it will cause before it actually causes it. For the pitch controlled wind turbine looked at here The desired signal to feed forward is the instant wind speed experienced by the rotor. This signal is very difficult to obtain on a standard wind turbine because the anemometer on top of the nacelle is not accurate enough due to the disturbance from the rotor itself and it is not practical to have an array of anemometers in all directions to get an upstream measurement of the wind speed.

[0150] The active control scheme can also be utilised for energy production optimisation. The pitch angle or rotor speed is then controlled in order to operate at maximum C_p. Various ways of implementing the optimisation are possible but the simplest is to use the wind speed and then, based on steady state characteristics, determine the optimal pitch angle or rotational speed. The principle for a variable speed machine can be illustrated as in FIG. 7. It is desired to operate the rotor at a constant tip speed ratio in order to have maximum efficiency. The tip speed ratio is the ratio between the wind speed experienced by the whole rotor and the rotational speed of the rotor. It is therefore very desirable to know the wind speed.

[0151] If the wind speed is to be used for control purposes there are several considerations that has to be taken into account. It is clear that the better the wind is known the safer it can be fed forward.

[0152] The wind speed that is needed for power limit control is the wind speed seen (or to be seen) by the whole rotor. This wind speed is very difficult to obtain, because the wind turbine itself disturbs the flow and the wind speed is measured at one point.

[0153] On the current wind turbines the cup anemometer is mounted on the nacelle. The wind speed measured there can not be used as wind speed input for power limit control because the wind speed at the anemometer is disturbed by the rotor itself. Instead of using an anemometer on top of the nacelle a laser anemometer can be used to measure the wind speed at a point (or several points) up stream.

[0154] The measured wind speed upstream can then be used to predict the wind speed that the rotor will see at some point later in time. The wind speed prediction will have to take the delay, the spatial filtering of the rotor and correlation between the wind speed measured at one point and the wind turbine's position. The estimate of the wind speed can be derived in two ways. It can be done in open loop as in FIG. 8 or it can be done in closed loop as in FIG. 9. The estimated wind speed can then be used in the control algorithm. All the models involved are dynamic. The wind turbine model should be as simple as possible only taking low frequency (˜1 Hz) dynamics into account.

[0155] Since the final goal of this investigation is to assess, with blade pitch control in mind, the possibility of using the upwind air velocity, measured from the wind turbine, a number of micrometeorological questions have to be answered. If we need five second of time, say, to carry out a change in the pitch, then we must know how well the measured wind speed, perhaps ten rotor diameters upstream, represents the actual wind speed at the rotor five seconds later.

[0156] In the following we will address the questions:

[0157] 1. Eddies of smaller extension than the linear dimensions of the rotor are uncorrelated over these dimensions and of no importance to the overall forcing. Consequently, we must emulate the area averaging of the wind field by the rotor as well as possible in the upwind detection. We assume that this is carried out by a line averaging in the wind direction. What should the length of this line then be?

[0158] 2. What is the correlation coefficient between the upstream wind speed and the wind speed at the rotor?

[0159] 3. How do we interpret the correlation coefficient to tell us how well the measured upstream wind speed matches the wind speed at the rotor.

[0160] Let the spectral tensor component for the streamwise, unfiltered velocity component u⁰(χ) be Φ_uu⁰(κ). Then spectral the tensor component of the spatially filtered velocity component u(χ) will be given by

[0161] Φ_uu(κ)=Ξ(κ)Φ_uu⁰(κ),  (4.1)

[0162] where the actual form of Ξ(κ) depends on the relation between u(χ) and u⁰(χ).

[0163] In general this relation is given by the integral over the entire physical space $\begin{array}{cc}u\left(x\right)=\underset{V}{\oint }{u}^0\left(x-x^{\prime }\right)h\left(x^{\prime }\right)d^3{x}^{\prime }.& \left(4.2\right)\end{array}$

[0164] The volume transfer function Ξ(κ) is then given by

[0165] Ξ(κ)=|ĥ(κ)|²,  (4.3)

[0166] where $\begin{array}{cc}\hat{h}\left(k\right)=\underset{V}{\oint }h\left(x\right)\exp \left(\mathrm{ik}·x\right)d^{3}x.& \left(4.4\right)\end{array}$

[0167] To obtain the usual one-dimensional spectrum F_u(κ₁) of u as a function of the streamwise component of the wave vector κ, we carry out the double integration over the other two components of κ, viz. $\begin{array}{cc}{F}_u\left(k_1\right)={\int }_{-\infty }^{\infty }d{k}_2{\int }_{-\infty }^{\infty }d{k}_3{\Phi }_{\mathrm{uu}}\left(k\right)={\int }_{-\infty }^{\infty }d{k}_2{\int }_{-\infty }^{\infty }d{k}_3\Xi \left(k\right){\Phi }_{\mathrm{uu}}^0\left(k\right).& \left(4.5\right)\end{array}$

[0168] The corresponding unfiltered one-dimensional spectrum F_u⁰(κ) is of course given by $\begin{array}{cc}{F}_u^0\left(k_1\right)={\int }_{-\infty }^{\infty }d{k}_2{\int }_{-\infty }^{\infty }d{k}_3{\Phi }_{\mathrm{uu}}^0\left(k\right).& \left(4.6\right)\end{array}$

[0169] The one-dimensional transfer function is defined as $\begin{array}{cc}H\left(k\right)=\frac{F_u\left(k\right)}{F_u^0\left(k\right)}.& \left(4.7\right)\end{array}$

[0170] In the special case of isotropy we have $\begin{array}{cc}{\Phi }_{\mathrm{uu}}^0\left(k\right)=\frac{E\left(k\right)}{4\pi }\frac{k_2^2+k_3^2}{k^4},& \left(4.8\right)\end{array}$

[0171] where E(κ) is the energy spectrum and in this case there is a simple relation between F_u⁰(κ) and E(κ), namely $\begin{array}{cc}E\left(k\right)=k^3\frac{d}{dk}\frac{1}{k}\frac{d{F}_u^0}{dk}.& \left(4.9\right)\end{array}$

[0172] In this case we assume that h(χ)=h_L(χ) is given by $\begin{array}{cc}{h}_L\left(x\right)=\frac{1}{\pi }\frac{L}{L^2+x_1^2}\delta \left(x_2\right)\delta \left(x_3\right),& \left(4.10\right)\end{array}$

[0173] where L, as indicated in FIG. 19, is the characteristic length of the one-dimensional spatial weighting. Substituting in (4.4), we get $\begin{array}{cc}{\hat{h}}_L\left(k\right)=e^{-\left|k_1\right|L}.& \left(4.11\right)\end{array}$

[0174] The volume transfer function thus becomes $\begin{array}{cc}{\Xi }_L\left(k\right)=e^{-2\left|k_1\right|L}.& \left(4.12\right)\end{array}$

[0175] Since Ξ_L(κ) is independent of κ₂and κ₃, (4.5), (4.6) and (4.7) imply that the one-dimensional transfer function is equal to H_L(κ), i.e. $\begin{array}{cc}{H}_L\left(k_1\right)={\Xi }_L\left(k\right)=e^{-2\left|k_1\right|L}.& \left(4.13\right)\end{array}$

[0176] We now assume that the mean wind and, consequently the χ₁-direction is perpendicular to the rotor of the wind turbine. The spatial weighting function h(χ)=h_R(χ) is therefore $\begin{array}{cc}{h}_R\left(x\right)=\delta \left(x_1\right)\times \{\begin{array}{cc}\frac{1}{\pi R^2},& \sqrt{x_2^2+x_3^2}<R\\ 0,& \sqrt{x_2^2+x_3^2}\ge R\end{array},& \left(4.14\right)\end{array}$

[0177] where R is the radius of the rotor (see FIG. 19).

[0178] We substitute (4.14) in (4.4) which because of the delta function reduces to the double integral $\begin{array}{cc}{\hat{h}}_R\left(k\right)=\underset{\mathrm{Rotor}}{\oint }\int d{x}_{2}d{x}_3\frac{\exp \left(i\left\{k_2{x}_2+k_3{x}_3\right\}\right)}{\pi R^2}.& \left(4.15\right)\end{array}$

[0179] This integral can easily be carried out by using plane polar integration variables: $\begin{array}{cc}\begin{array}{c}{\hat{h}}_R\left(k\right)=\frac{1}{\pi R^2}{\int }_0^{R}rdr{\int }_0^{2\pi }\exp \left(ir\left\{k_2\cos \varphi +k_3\sin \varphi \right\}\right)d\varphi \\ =\frac{1}{\pi R^2}{\int }_0^{R}rdr{\int }_0^{2\pi }\exp \left(ir\sqrt{k_2^2+k_3^2}\cos {\varphi }^{\prime }\right)d{\varphi }^{\prime }\\ =\frac{2}{R^2}{\int }_0^R{J}_0\left(\sqrt{k_2^2+k_3^2}r\right)rdr\\ =2\frac{J_1\left(\sqrt{k_2^2+k_3^2}R\right)}{\sqrt{k_2^2+k_3^2}R},\end{array}& \left(4.16\right)\end{array}$

[0180] where J_v(s) is the Bessel Function of the first kind and of order v.

[0181] The volume transfer function Ξ_R(κ) is according to (4.3) $\begin{array}{cc}{\Xi }_R\left(k\right)=4\frac{J_1^2\left(\sqrt{k_2^2+k_3^2}R\right)}{\left(k_2^2+k_3^2\right)R^2}.& \left(4.17\right)\end{array}$

[0182] Substituting this into (4.5), we obtain the filtered one-dimensional spectrum F_u^R(κ₁). In order to be able to carry out the two integrations we must assume a particular form of the spectral tensor component Φ_uu⁰(κ) and here we assume that the turbulent velocity field is isotropic so that we can use (4.8). The integration could be carried out over κ₂and κ₃, but because of the form of H_R(κ) in (4.17) and χ_uu^{0(κ) in (}4.8) it is

[0183] advantageous to change to plane polar integration variables $\begin{array}{cc}\left\{\begin{array}{c}{k}_2\\ k_3\end{array}\right\}=K\times \left\{\begin{array}{c}\cos \theta \\ \sin \theta \end{array}\right\}.& \left(4.18\right)\\ \begin{array}{c}{F}_u^R\left(k_1\right)={\int }_{-\infty }^{\infty }d{k}_2{\int }_{-\infty }^{\infty }d{k}_3{\Xi }_R\left(k\right){\Phi }_{\mathrm{uu}}^0\left(k\right)\\ ={\int }_{-\infty }^{\infty }d{k}_2{\int }_{-\infty }^{\infty }d{k}_{3}4\frac{J_1^2\left(\sqrt{k_2^2+k_3^2}R\right)}{\left(k_2^2+k_3^2\right)R^2}\frac{E\sqrt{k_1^2+k_2^2+k_3^2}}{4\pi }{k}_2^2+\frac{k_3^2}{{\left(k_1^2+k_2^2+k_3^2\right)}^2}\\ =\frac{2}{R^2}{\int }_0^{\infty }{J}_1^2\left(\mathrm{KR}\right)\frac{E\left(\sqrt{k_1^2+K^2}\right)}{{\left(k_1^2+K^2\right)}^2}KdK.\end{array}& \left(4.19\right)\end{array}$

[0184] We cannot carry out the final integration over K without making a specific assumption about the energy spectrum E(κ) and here we chose the simplest possible, namely that

[0185] E(κ)=α∈^2/3κ^−5/3, 0<κ<∞,  (4.20)

[0186] where α≈1.7 is the dimensionless Kolmogorov constant and ∈ the rate of dissipation of specific kinetic energy.

[0187] Equation (4.9) implies that $\begin{array}{cc}{F}_u^0\left(k\right)=\frac{9}{55}E\left(\left|k\right|\right)=\frac{9}{55}\alpha {ϵ}^{2/3}{\left|k\right|}^{-5/3}.& \left(4.21\right)\end{array}$

[0188] Substituting (4.20) in (4.19), we get $\begin{array}{cc}\begin{array}{c}{F}_u^R\left(k\right)=\frac{2}{R^2}{\int }_0^{\infty }\frac{{\mathrm{KJ}}_1^2\left(\mathrm{KR}\right)}{{\left(k^2+K^2\right)}^2}\alpha {{ϵ}^{2/3}\left(k^2+K^2\right)}^{-5/6}dK\\ =\alpha {ϵ}^{2/3}{R}^{5/3}2{\int }_0^{\infty }\frac{{\mathrm{sJ}}_1^2\left(s\right)}{{\left(k^2{R}^2+s^2\right)}^{17/6}}ds\\ =\frac{9}{55}{\mathrm{\alpha ϵ}}^{2/3}{\left|k\right|}^{-5/3}\times H_R\left(k\right),\end{array}\mathrm{where}& \left(4.22\right)\\ H_R\left(k\right)\equiv \frac{F_u^R\left(k\right)}{F_u^0\left(k\right)}={}_{1}F_2\left(\frac{3}{2};\frac{1}{6},3;{\left(\mathrm{kR}\right)}^2\right)-{\left(\mathrm{kR}\right)}^{5/3}\frac{4}{\sqrt{\pi }}\frac{\Gamma \left(1/6\right)\Gamma \left(7/3\right)}{\Gamma \left(11/6\right)\Gamma \left(23/6\right)}{}_{1}F_2\left(\frac{7}{3};\frac{11}{6},\frac{23}{6};{\left(\mathrm{kR}\right)}^2\right),& \left(4.23\right)\end{array}$

[0189] where ₁F₂(a;b,c;χ) is the Generalized Hypergeometric Function.

[0190] In the derivation we have assumed uniform weighting over the entire rotor plane. However, since the largest part of the external torque from the wind is located far from the axis, we have also calculated H_R(κ) in the case when the spatial weighting over the rotor plane increases with the square of the distance from the rotor axis. Again, it was possible to derive an explicit expression like (4.23), albeit a lot more complicated.

[0191] FIG. 20 shows H_R(κ), given by (4.23) together with a simple exponential expression exp(−0.75κR). Also shown is H_R(κ) in the case when the uniform weighting over the rotor is replaced with a weighting which increased with the square of the distance from the rotor axis.

[0192] Assuming that the exponential expression is accurate enough to characterize the area filtering by the rotor, we conclude that the lidar line averaging should have a characteristic length L of about

[0193] L≈0.375R,  (4.24)

[0194] irrespective of the weighting over the rotor plane.

[0195] The question about how well correlated the upstream longitudinal velocity component is to the same velocity component just in front of the rotor of the wind turbine will be discussed in terms of the so-called longitudinal velocity coherence. If Taylor's hypothesis were strictly valid this coherence would be exactly one and we would expect perfect correlation. However, Taylor's hypothesis is of course only approximately true and we must consequently move a little outside the realm of standard assumptions when we deal with the longitudinal coherence.

[0196] From an experimental point of view, deviations from unity longitudinal coherence are difficult to determine. The reason is that, even when making a strong effort, the line between the two observations points will in general only by coincidence be sufficiently close to mean-wind direction during the observation period. This direction is known only after the measurements have been carried out and if there is just a small angle between the mean-wind vector and the separation line between the anemometers the measured coherence will be dominated by lateral coherence which falls off with frequency much faster than the longitudinal coherence. We will therefore first make an attempt to understand the longitudinal coherence and the corresponding correlation on basis of a simple, semi-quantitative model and then estimate the contamination from lateral coherence due to misalignment.

[0197] Consider the covariance

[0198] R_u(X,τ)=((u(χ,t)−U)u(χ+X, t+τ))−U)),  (4.25)

[0199] where u(χ,t) is the instantaneous wind-velocity component along the mean-wind direction at point χ and time t, U=(u(χ,t)) the mean-wind speed, X the separation vector and τ the time lag. Then the cross-spectrum is $\begin{array}{cc}{\chi }_u\left(X,\omega \right)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{R}_u\left(X,\tau \right)e^{-i\omega \tau }d\tau & \left(4.26\right)\end{array}$

[0200] and the coherence defined as $\begin{array}{cc}{\mathrm{coh}}_u\left(X,\omega \right)={\left|\frac{{\chi }_u\left(X,\omega \right)}{{\chi }_u\left(0,\omega \right)}\right|}^2.& \left(4.27\right)\end{array}$

[0201] This definition is general and can be used irrespective of the orientation of the separation vector X.

[0202] In the case when X is along the mean-wind direction, i.e. X=Xi₁, Taylor's hypothesis implies $\begin{array}{cc}{R}_u\left(X{i}_1,\tau \right)=R_u\left(0,\tau -\frac{X}{U}\right)& \left(4.28\right)\end{array}$

[0203] and, when substituting this in (4.26), we obtain $\begin{array}{cc}{\chi }_u\left(X{i}_1,\omega \right)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{R}_u\left(0,\tau -\frac{X}{U}\right)e^{-i\omega \tau }d\tau =e^{-i\frac{X}{U}\omega }{\chi }_u\left(0,\omega \right).& \left(4.29\right)\end{array}$

[0204] We see that there is only a phase difference between χ_u(Xi₁,ω) and χ_u(0,ω). Consequently the coherence, defined by (4.27), must be one, irrespective of the frequency, in cases where we assume Taylor's hypothesis.

[0205] Another way of stating this consequence of Taylor's hypothesis is to interpret (4.28) and note that u((χ+X)i₁, t) is perfectly correlated with u(χi₁,t−X/U)= the upstream, streamwise velocity component at a time so much earlier as it takes a disturbance to travel the distance X with the speed U. In other words, $\begin{array}{cc}\begin{array}{c}⟨\left\{u\left(x{i}_1,t-\frac{X}{U}\right)-U\right\}\left\{u\left(\left[x+X\right]i_1,t\right)-U\right\}⟩=R_u\left(X{i}_1,\frac{X}{U}\right)\\ =R_u\left(0,0\right)={\sigma }_u^2,\end{array}& \left(4.30\right)\end{array}$

[0206] where σ_u²is the variance of u.

[0207] We will modify (4.29) to account for the decorrelation due to the turbulence by the generalization (Kristensen, 1979) $\begin{array}{cc}{\chi }_u\left(X{i}_1,\omega \right)=\exp \left(-\frac{X}{\mathrm{UT}\left(\omega \right)}\right)e^{-i\frac{X}{U}\omega }{\chi }_u\left(0,\omega \right),& \left(4.31\right)\end{array}$

[0208] where the exponential ‘decay factor’ is a decreasing function of X. The ‘eddy mean-lifetime’τ(ω) is assumed to be a function of the eddy size λ˜2πU/ω^↑.

[0209] ↑Actually we cheat a little here and in the following and use Taylor's hypothesis to translate eddy size or reciprocal wave number to frequency ω.

[0210] Mann (1994) found that $\begin{array}{cc}T\left(\omega \right)=\Gamma \kappa \frac{z}{u_{*}}{\left(\frac{\left|\omega \right|}{U}z\right)}^{-2/3}{\left\{{}_{2}F_1\left(\frac{1}{3},\frac{17}{6};\frac{4}{3};-{\left[\frac{\left|\omega \right|}{U}z\right]}^{-2}\right)\right\}}^{-1/2},& \left(4.32\right)\end{array}$

[0211] where Γ≈3 is a dimensionless constant, κ≈0.4 the von Karman constant, u. the friction velocity and z the height, brings his spectral model in good agreement with experimental data over land and over water.

[0212] for very numerically very small and very large values of ω/U we have $\begin{array}{cc}T\left(\omega \right)=\Gamma \kappa \frac{z}{u_{*}}\times \{\begin{array}{cc}\sqrt{\frac{2}{3B\left(\frac{3}{2},\frac{4}{3}\right)}}{\left(\frac{\left|\omega \right|}{U}z\right)}^{-1}\approx 1.2{\left(\frac{\left|\omega \right|}{U}z\right)}^{-1},& \left|\omega \right|\to \infty \\ {\left(\frac{\left|\omega \right|}{U}z\right)}^{-2/3},& \left|\omega \right|\to 0,\end{array}\mathrm{where}& \left(4.33\right)\\ B\left(x,y\right)=\frac{\Gamma \left(x\right)\times \Gamma \left(y\right)}{\Gamma \left(x+y\right)}& \left(4.34\right)\end{array}$

[0213] is the (Complete) Beta Function.

[0214] As we are concerned with rather large eddies, we shall use the asymptote $\begin{array}{cc}T\left(\omega \right)\approx 1.5\frac{z}{u_{*}}{\left(\frac{\omega }{U}z\right)}^{-1}.& \left(4.35\right)\end{array}$

[0215] All we now have to do to determine the correlation coefficient with longitudinal separation $\begin{array}{cc}\rho \left(X{i}_1\right)=\frac{R_u\left(X{i}_1,\frac{X}{U}\right)}{{\sigma }_u^2}& \left(4.36\right)\end{array}$

[0216] is to apply the transform $\begin{array}{cc}{R}_u\left(X{i}_1,\tau \right)={\int }_{-\infty }^{\infty }{\chi }_u\left(X{i}_1,\omega \right)e^{i\omega \tau }d\omega ,& \left(4.37\right)\end{array}$

[0217] with τ=X/U, to χ_u(Xi₁,ω) given by (4.31): $\begin{array}{cc}\begin{array}{c}{R}_u\left(X{i}_1,\frac{X}{U}\right)={\int }_{-\infty }^{\infty }\exp \left(-\frac{2}{3}\frac{u_{*}}{U}\frac{X}{z}\left\{\frac{\omega }{U}z\right\}\right){\chi }_u\left(0,\omega \right)d\omega \\ ={\int }_{-\infty }^{\infty }\exp \left(-\frac{2}{3}\frac{u_{*}}{U}\frac{X}{z}\mathrm{kz}\right)F_u\left(k\right)dk,\end{array}& \left(4.38\right)\end{array}$

[0218] where F_u(κ) is defined by (4.5).

[0219] We assume, according (see FIG. 20) that

[0220] F_u(κ)=e^−0.75κRF_u⁰(κ)  (4.39)

[0221] and we will use the von Kármán spectrum with a length scale equal to z, consistent with Mann (1994), i.e. $\begin{array}{cc}{F}_u^0\left(k\right)=\frac{{\sigma }_u^2}{B\left(1/2,1/3\right)}\frac{z}{{\left(1+{\left(\mathrm{kz}\right)}^2\right)}^{5/6}}.& \left(4.40\right)\end{array}$

[0222] Further, we will assume that the rotor radius R is half the hub height z.

[0223] Substituting in (4.38), we get $\begin{array}{cc}{R}_u\left(X{i}_1,\frac{X}{U}\right)=\frac{2{\sigma }_u^2}{B\left(1/2,1/3\right)}G\left(\frac{1}{3}\frac{u_{*}}{U}\frac{X}{R}+\frac{3}{8}\right),\mathrm{where}& \left(4.41\right)\\ \begin{array}{c}G\left(\xi \right)={\int }_0^{\infty }\frac{e^{-\xi s}}{{\left(1+s^2\right)}^{5/6}}ds\\ =3{\xi }_1{F}_2\left(1;\frac{7}{6},\frac{3}{2};-{\left\{\frac{\xi }{2}\right\}}^2\right)-\frac{{\pi }^{3/2}}{\Gamma \left(5/6\right)}{\left(\frac{\xi }{2}\right)}^{1/3}{Y}_{-1/3}\left(\xi \right),\end{array}& \left(4.42\right)\end{array}$

[0224] where Y_v(s) is the Bessel Function of the second kind and of order v. is We can now determine the correlation coefficient ρ(Xi₁): $\begin{array}{cc}\rho \left(X{i}_1\right)=\frac{R_u\left(X{i}_1,\frac{X}{U}\right)}{R_u\left(0,0\right)}=\frac{G\left(\frac{1}{3}\frac{u_{*}}{U}\frac{X}{R}+\frac{3}{8}\right)}{G\left(\frac{3}{8}\right)}.& \left(4.43\right)\end{array}$

[0225] The correlation function is shown in FIG. 21.

[0226] If the normal to the rotor plane, the direction of the yaw, is not parallel to the mean wind direction, a gust detected at the distance X will not travel along the yaw direction and will therefore not hit the center of the rotor. The correlation will consequently not be perfect between the upstream wind velocity and the velocity at the rotor at the time X/U later. Let the angle between the yaw direction and the mean-wind direction be α. Then, assuming Taylor's hypothesis, the correlation in question is that between velocity components laterally separated the distance X tan(α)≈Xα. By simple modelling we can estimate this correlation and we will consider both point velocities and velocities which are averaged over the rotor plane. We start with point velocities.

[0227] The spectral tensor component for the longitudinal velocity is given by (4.8) and here we use the energy spectrum which can be derived from the von Kármán spectrum and the relation (4.9): $\begin{array}{cc}E\left(k\right)=\frac{2}{3}\frac{{\sigma }_u^{2}z}{B\left(3/2,4/3\right)}\frac{{\left(\mathrm{kz}\right)}^4}{{\left(1+{\left(\mathrm{kz}\right)}^2\right)}^{17/6}}.& \left(4.44\right)\end{array}$

[0228] The spectral tensor component can be written $\begin{array}{cc}{\Phi }_{\mathrm{uu}}^0\left(k\right)=\frac{1}{6\pi }\frac{{\sigma }_u^2{z}^3}{B\left(3/2,4/3\right)}\frac{\left(k_2^2+k_3^2\right)z^2}{{\left\{1+\left(k_1^2+k_2^2+k_3^2\right)z^2\right\}}^{17/6}}.& \left(4.45\right)\end{array}$

[0229] The covariance of the two unfiltered velocity components with the lateral separation Xα becomes $\begin{array}{cc}\begin{array}{c}{R}_u^0\left(X\alpha i_2\right)=∯\int {\Phi }_{\mathrm{uu}}^0\left(k\right)e^{i{k}_{2}X\alpha }{d}^{3}k\\ =\frac{1}{6\pi }\frac{{\sigma }_u^2{z}^3}{B\left(3/2,4/3\right)}{\int }_{-\infty }^{\infty }d{k}_1{\int }_0^{\infty }KdK\\ {\int }_0^{2\pi }d\phi K^2{z}^2\frac{\exp \left(i\mathrm{KX}\alpha \cos \phi \right)}{{\left\{1+\left(k_1^2+K^2\right)z^2\right\}}^{17/6}}\\ =\frac{2}{3}\frac{{\sigma }_u^2}{B\left(3/2,4/3\right)}{\int }_0^{\infty }\frac{dt}{{\left(1+s^2+t^2\right)}^{17/6}}{\int }_0^{\infty }{J}_0\left(\frac{X\alpha }{z}s\right)s^{3}ds\\ =\frac{8}{9}{\sigma }_u^2{\int }_0^{\infty }{J}_0\left(\frac{X\alpha }{z}s\right)\frac{s^3}{{\left(1+s^2\right)}^{7/3}}ds={\sigma }_u^2{\rho }^0\left(X\alpha i_2\right),\end{array}& \left(4.46\right)\end{array}$

[0230] where ρ⁰(Xαi₂) is the correlation function, given by $\begin{array}{cc}{\rho }^0\left(X\alpha i_2\right)=\frac{8}{9}{\sigma }_u^2\left\{\frac{1}{\Gamma \left(4/3\right)}{\left(\frac{X\alpha }{2z}\right)}^{1/3}{K}_{1/3}\left(\frac{X\alpha }{z}\right)-\frac{1}{\Gamma \left(7/3\right)}{\left(\frac{X\alpha }{2z}\right)}^{4/3}{K}_{4/3}\left(\frac{X\alpha }{z}\right)\right\},& \left(4.47\right)\end{array}$

[0231] where K_v(s) is the a Modified Bessel Function (convergent for s→∞) of order v. Here and in the following we consider the typical situation where the hub height is twice the radius of the rotor, i.e. z=2R. The last result can therefore be written $\begin{array}{cc}{\rho }^0\left(X\alpha i_2\right)=\frac{8}{9}{\sigma }_u^2\left\{\frac{1}{\Gamma \left(4/3\right)}{\left(\frac{X\alpha }{4R}\right)}^{1/3}{K}_{1/3}\left(\frac{X\alpha }{2R}\right)-\frac{1}{\Gamma \left(7/3\right)}{\left(\frac{X\alpha }{4R}\right)}^{4/3}{K}_{4/3}\left(\frac{X\alpha }{2R}\right)\right\}.& \left(4.48\right)\end{array}$

[0232] When the velocity signal is averaged over the rotor we must replace Φ_uu⁰(κ) by $\begin{array}{cc}{\Phi }_{\mathrm{uu}}\left(k\right)=2^5\frac{2}{3\pi }\frac{{\sigma }_u^2{R}^3}{B\left(3/2,4/3\right)}\frac{J_1^2\left(\sqrt{k_2^2+k_3^2}R\right)}{{\left\{1+4\left(k_1^2+k_2^2+k_3^2\right)R^2\right\}}^{17/6}}.& \left(4.49\right)\end{array}$

[0233] The covariance becomes $\begin{array}{cc}\begin{array}{c}{R}_u\left(X\alpha i_2\right)=∯\int {\Phi }_{\mathrm{uu}}\left(k\right)e^{i{k}_{2}X\alpha }{d}^{3}k\\ =2^5\frac{2}{3\pi }\frac{{\sigma }_u^2{R}^3}{B\left(3/2,4/3\right)}{\int }_{-\infty }^{\infty }d{k}_1{\int }_0^{\infty }KdK\\ {\int }_0^{2\pi }d\phi J_1^2\left(\mathrm{KX}\alpha \right)\frac{\exp \left(i\mathrm{KX}\alpha \cos \phi \right)}{{\left\{1+4\left(k_1^2+K^2\right)R^2\right\}}^{17/6}}\\ =2^6\frac{4}{9}{\sigma }_u^2{\int }_0^{\infty }\frac{J_1^2\left(s\right)J_0\left(\left[X\alpha /R\right]s\right)}{{\left\{1+4s^2\right\}}^{7/3}}sds.\end{array}& \left(4.50\right)\end{array}$

[0234] The correlation coefficient $\begin{array}{cc}\rho \left(X\alpha i_2\right)={\int }_0^{\infty }\frac{J_1^2\left(s\right)J_0\left(\left[X\alpha /R\right]s\right)}{{\left\{1+4s^2\right\}}^{7/3}}sds/{\int }_0^{\infty }\frac{J_1^2\left(s\right)}{{\left\{1+4s^2\right\}}^{7/3}}sds& \left(4.51\right)\end{array}$

[0235] can be evaluated by numerical integration as a function of Xα/R. FIG. 22 shows the correlation when the velocity is unfiltered and when it is area filtered.

[0236] Until now we have only discussed the correlation coefficient ρ₀between the upstream velocity and the velocity at the wind turbine. It seems tempting to assume that if ρ₀is close to one then a large turbulent excursion of the upstream velocity will give rise to an excursion which is similar in shape and magnitude at the wind turbine at a time X/U later. How can this assumption be quantified?

[0237] Let us call the upstream velocity u₁and the velocity at the turbine u. Let the variance σ_u²and the correlation coefficient ρ₀be given, and let us assume that u₁and u have a joint Gaussian probability density $\begin{array}{cc}P\left(u_1,u\right)=\frac{\exp \left(-\frac{1}{2\left(1-{\rho }_0^2\right){\sigma }_u^2}\left\{u_1^2-2{\rho }_0{u}_{1}u+u^2\right\}\right)}{2\pi \sqrt{1-{\rho }_0^2}{\sigma }_u^2}.& \left(4.52\right)\end{array}$

[0238] We may recast this expression as a product of the probability density function ρ₀(u₁) for u₁and the conditional probability density P₁(u₁;u) for u. We are particular interested in P₁which becomes $\begin{array}{cc}{P}_1\left(u_1;u\right)=\frac{1}{\sqrt{2\pi }\sqrt{1-{\rho }_0^2}{\sigma }_u}\exp \left(-\frac{{\left(u-{\rho }_0{u}_1\right)}^2}{2\left(1-{\rho }_0^2\right){\sigma }_u^2}\right).& \left(4.53\right)\end{array}$

[0239] We see from this expression that the conditional mean value of u is ρ₀u₁and not u₁. The implication is that if ρ₀<1, the mean (and the most probable) value of the velocity u at the wind turbine is less than u₁if u₁is positive. A given upstream excursion will be suppressed on its way to the wind turbine.

[0240] Since u is biased in such a way that its most probable value is u₀=ρ₀u₁, it seems reasonable to use this value for controlling the rotor. Consequently, we may ask the question: what is the probability ψ₁that |u−u₀|/|u₀| is less than a small positive number ∈?

[0241] This probability will depend on ρ₀, ∈ and the excursion in units of the standard deviation u₀/σ_u. $\begin{array}{cc}\begin{array}{c}{\psi }_1\left({\rho }_0,\epsilon ,\frac{u_0}{{\sigma }_u}\right)={\int }_{u_0-\epsilon \left|u_0\right|}^{u_0+\epsilon \left|u_0\right|}{P}_1\left(u_1,u\right)du\\ =\mathrm{erf}\left(\frac{\epsilon }{\sqrt{2\left(1-{\rho }_0^2\right)}}\frac{\left|u_0\right|}{{\sigma }_u}\right),\end{array}\mathrm{where}& \left(4.54\right)\\ \mathrm{erf}\left(x\right)=\frac{2}{\sqrt{\pi }}{\int }_0^x{e}^{-s^2}ds& \left(4.55\right)\end{array}$

[0242] is the Error Function, shown in FIG. 23.

[0243] When χ>2, the Error Function can conveniently be approximated by $\begin{array}{cc}\mathrm{erf}\left(x\right)\approx 1-\frac{e^{-x^2}}{\sqrt{\pi }x}.& \left(4.56\right)\end{array}$

[0244] We might get a better estimate of the excursion u at the wind turbine if we measured the wind speed at several places upstream (a range-gating system). Let us assume that we measure the wind speeds at two upstream positions and detect the excursions u₁and u₂. Let the correlation coefficients between u and u₁and u₂be ρ₁and ρ₂, respectively, and let the correlation coefficient between u₁and u₂be ρ. Further, let us assume that u, u₁and u₂have a joint Gaussian probability density. Then the conditional probability for u, when u₁and u₂are given, becomes $\begin{array}{cc}{P}_2\left(u_1,u_2;u\right)=\frac{1}{\sqrt{2\pi }}\sqrt{\frac{1-{\rho }^2}{1+2{\mathrm{\rho \rho }}_1{\rho }_2-{\rho }^2-{\rho }_1^2-{\rho }_2^2}}\frac{1}{{\sigma }_u}\times \exp \left(-\frac{1-{\rho }^2}{1+2{\mathrm{\rho \rho }}_1{\rho }_2-{\rho }^2-{\rho }_1^2-{\rho }_2^2}\frac{1}{2{\sigma }_u^2}{\left\{u-\left[\frac{{\rho }_1-{\mathrm{\rho \rho }}_2}{1-{\rho }^2}{u}_1+\frac{{\rho }_2-{\mathrm{\rho \rho }}_1}{1-{\rho }^2}{u}_2\right]\right\}}^2\right).& \left(4.57\right)\end{array}$

[0245] We see that in this case the most probable value of u becomes $\begin{array}{cc}{u}_0=\frac{{\rho }_1-{\mathrm{\rho \rho }}_2}{1-{\rho }^2}{u}_1+\frac{{\rho }_2-{\mathrm{\rho \rho }}_1}{1-{\rho }^2}{u}_2.& \left(4.58\right)\end{array}$

[0246] The probability ψ₂that |u−u₀| is less than ∈|u₀|is in this case given by $\begin{array}{cc}{\psi }_2\left({\rho }_1,{\rho }_2,\rho ,\epsilon ,\frac{u_0}{{\sigma }_u}\right)=\mathrm{erf}\left(\sqrt{\frac{1-{\rho }^2}{2\left(1+2{\mathrm{\rho \rho }}_1{\rho }_2-{\rho }^2-{\rho }_1^2-{\rho }_2^2\right)}}\epsilon \frac{\left|u_0\right|}{{\sigma }_u}\right).& \left(4.59\right)\end{array}$

[0247] In the limit when ρ=1 and ρ₁=ρ₂≡ρ₀we obtain the same result as in the case when there is only one upstream measurement with u₀=ρ₀(u₁+u₂)/2.

[0248] It might be interesting to see if we get a better estimate with two upstream measurements if ρ₁=ρ₂≡ρ₀and ρ<1. In this case (4.59) becomes $\begin{array}{cc}{\psi }_2\left({\rho }_0,{\rho }_0,\rho ,\epsilon ,\frac{u_0}{{\sigma }_u}\right)=\mathrm{erf}\left(\sqrt{\frac{1+\rho }{2\left(1+\rho -2{\rho }_0^2\right)}}\epsilon \frac{\left|u_0\right|}{{\sigma }_u}\right).\mathrm{Since}& \left(4.60\right)\\ \frac{1+\rho }{2\left(1+\rho -2{\rho }_0^2\right)}\ge \frac{1}{2\left(1-{\rho }_0^2\right)}& \left(4.61\right)\end{array}$