- Email: [email protected]
Simulation of aggregate wind farm short-term production variations
Renewable Energy 35 (2010) 2602e2609
Contents lists available at ScienceDirect
Renewable Energy journal homepage: www.elsevier.com/locate/renene
Simulation of aggregate wind farm short-term production variations R. Goi c, J. Krstulovi c*, D. Jakus Faculty of Electrical Engineering and Naval Architecture, University of Split, Rudera Boskovica bb, 21000 Split, Croatia
a r t i c l e i n f o
a b s t r a c t
Article history: Received 25 November 2009 Accepted 7 April 2010 Available online 1 May 2010
The variability of wind power production poses the greatest challenge in the integration of large-scale wind power in power systems. Furthermore, larger-scale penetration implies a wider geographical spreading of installed wind power, resulting in reduced variability and the smoothing effect of total power generation. Therefore, analysis of the impact of wind power variations on power system operation requires adequate modeling of aggregate power output from geographically dispersed wind farms. This paper analyzes different aspects of Markov chain Monte Carlo simulation methods for the synthetic generation of dependent wind power time series. However, testing indicates that these approaches do not adequately model the stochastic dependence between wind power time series in conjunction with individual persistence which is necessary to obtain realistic distributions of aggregate power output and total power variations. Consequently, a novel approach based on a modified second-order Markov chain Monte Carlo simulation is proposed. Simulation results show that this method obtains synthetic time series of aggregate wind power which very closely fit the original data, with respect to both the cumulative density function of total output power and the probability density function of power variations. Ó 2010 Elsevier Ltd. All rights reserved.
Keywords: Aggregate wind power generation Wind power variations Monte Carlo simulation Markov chain Stochastic dependency
1. Introduction With an installed nominal capacity of 121.2 GW, wind power makes about 1.5% of world wide electricity generation (up to the end of 2008). This number is growing rapidly, having doubled in the three years between 2005 and 2008. Several countries have achieved relatively high levels of wind penetration, such as 19% of electricity production in Denmark, 11% in Spain and Portugal, and 7% in Germany in 2008 [1]. With this trend of growing wind power capacities, the integration of large-scale wind power in power systems has become more difficult. The stochastic nature of wind creates fluctuating and intermitting wind power. Since instantaneous electrical generation and consumption must remain balanced to maintain grid stability, such fluctuations create substantial obstacles in incorporating large amounts of wind power into a grid. Variations of wind power production occur on all timescales, from seconds to years, where even short-term variations are, to some extent, unpredictable. Observing the power variations at a single wind turbine or a single wind farm may lead to false conclusions. Namely, the geographical spreading of wind power reduces variability, increases predictability and reduces occasions with zero or peak output.
* Corresponding author. Tel.: þ385 21 305 802. E-mail address: [email protected] (J. Krstulovi c). 0960-1481/$ e see front matter Ó 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.renene.2010.04.005
Consequently, as wind farms are spread over wider geographical areas working under different wind conditions, the resulting total wind power generation becomes smoother. Designing a wind farm requires installing a meteorological station at a potential location and recording all the relevant wind data series. These time series are usually not long enough to create a representative model for potential wind farm generation and may lead to biased conclusions regarding energy capacity estimation, cost effectiveness, system design issues, power variations, etc. Consequently, it is necessary to use stochastic simulation techniques to synthetically generate wind speed (wind power) data while preserving its statistical parameters. Since wind speed is conditioned with lots of climatologic, atmospheric and environmental parameters, such phenomena are impossible to describe sufficiently using a deterministic approach. Furthermore, when studying the effects of varying wind power production on the power system, it is important to consider larger areas, i.e. aggregate power output of multiple wind farms. This additionally complicates the problem since simulating wind power for two or more sites also requires modeling their stochastic dependence in order to obtain realistic total power generation data and its variations. In related literature, various probabilistic and stochastic methods are used for synthetic generation of time series. The most relevant methods include: Monte Carlo simulation of independent random numbers distributed according to a Weibull probability distribution function, autoregressive (AR) in combination with
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
moving average (ARMA) models, Markov chain Monte Carlo (MCMC) methods, and a Wavelet e based approach. In the case when a certain persistency is required in a generated time series, the most commonly used methods are ARMA and MCMC. The basic advantage of the AR and ARMA models is controlled synthetization of a desired autocorrelation function (acf) while defining a relatively small number of parameters. On the other hand, the probability density function (pdf) of a time series generated with ARMA, generally does not match that of the measured data, leading to a wrong estimation of the energy yield [2]. Appropriate use of the MCMC method for wind speed (power) time series simulation leads to satisfactory pdf and acf matching between the generated and measured data. In this method, wind activity is modeled as a Markov chain, i.e. a discrete stochastic process with defined probabilities for transitions between states. In recent works [2e5], different aspects of using the MCMC method for the simulation of wind speed data are analyzed. G. Papaefthymiou and B. Klöckl [2] show that obtaining a stochastic model directly in the wind power domain is superior to indirect modeling of wind power using wind speeds. This leads to a reduced number of states requiring a lower order Markov chain at equal power data resolution. In this paper, emphasis is given to the simulation of aggregate wind farm power output with correct representation of their total power and its short term (10 min) variations. Wind power time series are independently generated for two potential wind farm sites using the MCMC simulation method. Analysis of the aggregate power output distribution and short-term variations indicates that certain dependency between the generated time series is required to obtain satisfactory results. Furthermore, to describe their dependency, an empirical copula is implemented in the MCMC simulation, resulting in good representation of aggregate wind power, but lacking the variations. To overcome these drawbacks, an algorithm is proposed using a modified MCMC simulation for synthetic generation of aggregate wind farm power output. Simulation results show that this method provides realistic representation of both aggregate power and its short term variations distribution. All the simulation methods were tested using 10 min average wind speed data over one year from two potential wind farm sites in Croatia. 2. Wind power variations The stochastic nature of wind causes significantly variable and intermittent wind turbine power output. This variability is present on all timescales: seconds, minutes, hours, days, months and years. Successfully defining and forecasting such variations makes crucial step for wind power integration and achieving optimal operation of a power system. The variability of the horizontal winds, that wind turbines interact with, shows interesting general characteristics approved with their kinetic energy of the speed fluctuation distribution in the frequency domain, known as the van der Hoven spectrum. Van der Hoven [6] shows that independently of the site, the typical spectrum exhibits two peaks approximately at 0.01 cycles/ hour (4 days period), and 50 cycles/h (1 min period), which are separated by an energy gap between periods of 10 min and 2 h. In a 4 days period, the peak corresponds to wind speed fluctuations due to the migratory pressure systems of a synoptic weather-map scale. The other excerpted peak is in the micrometeorological range as a result of the turbulence associated to local winds. The reason for the spectral gap is in the lack of a physical process which could support wind speed fluctuations in that frequency range. When comparing wind fluctuations of 10 min and hourly
2603
averages, both are placed in the spectral gap, whereas energy is higher in the 10 min period, closer to the high frequency peak, while in the hourly period the spectral gap mostly reaches its minimum. In this paper, only short-term power variations are analyzed, which include timescales from minutes to hours. Second-tosecond variations of individual wind turbine power output are insignificant due to the inertia of large rotating blades and the varying speed of the rotor. If the total power generated by a single wind farm is considered as a whole, variations of several seconds become absorbed since the same gusts do not occur simultaneously at all turbines. Furthermore, the total wind power variations of multiple wind farms are not significant, even for a few minutes period, depending on their geographical spreading. Namely, in larger geographical areas, weather fronts can not pass through all wind farms simultaneously. Thus, longer time scales are subject to the so-called smoothing effect, as illustrated in Fig. 1. Consequently, the maximum amplitudes of wind power variations experienced in the power system are reduced. Fig. 2 shows that higher variations of 10 min wind power occur less frequently when the total power of three wind farms is considered, in comparison to that of a single wind farm. This effect would be much more significant for a higher number of wind farms over an even larger geographic area. For better understanding of the geographical dispersion effect on short term wind power variations, it is important to observe the linear correlation coefficient between two wind power time series. Namely, dependence is most often measured with a linear correlation coefficient. This coefficient generally indicates the strength and direction of a linear relationship between two random variables and is calculated by:
P
ðx xÞðy yÞ
ffi rx;y ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi P P ðx xÞ2
ðy yÞ2
(1)
The linear correlation coefficient can range from 1 to 1, where the two extreme cases and the middle value are the most interesting to observe. For a maximal value of the correlation coefficient r ¼ 1, increases and decreases in power production from individual outputs occur simultaneously, resulting in the largest possible variations of total power production. For the case where r ¼ 0, the power outputs fluctuate independently, causing smaller variations. For the third case, when the correlation coefficient is r ¼ 1, the outputs fluctuate in opposite directions, resulting in minimal total
Fig. 1. An example of the smoothing effect by geographic dispersion for three wind farm locations in Croatia.
2604
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
3. Markov chain Monte Carlo simulation Available data series recorded at real wind farm locations are generally not long enough to be representative and sufficient for a detailed analysis. Thus, various simulation techniques are used for synthetic generation of new data. In this paper, the Markov chain Monte Carlo simulation method is applied. Using this approach, wind speed data is modeled as a Markov chain, while the Monte Carlo simulation is used for generating new data. The main objective of this simulation is to obtain realistic wind power output time series and can be obtained using one of two approaches: e first building a model in the wind speed domain and then transforming generated wind speeds to wind power or, e first transforming available wind speed data into wind power and then directly generating wind power time series. Fig. 2. The probability density function of 10 min wind power variations e the effect of geographical dispersion.
power variations. Therefore, smaller variations are obtained for lower correlation coefficients. In reality, this coefficient ranges between zero and one, depending on the geographical distance between the observed sites with inverse exponential dependency. Furthermore, considering multiple wind farm sites, correlation between sites could be represented using a correlation square matrix whose size is equal to number of sites. Each row of the correlation matrix represents the correlation between one site and all the other sites. Thus, the diagonal is composed of all ones while the remaining elements contain all cross-correlations between sites. Generally, fast power variations (within a minute time interval) are usually very low, with no impact on power system operation. On the other hand, in-hour variations (10e60 min) can be much more significant, with extreme ramp rates reaching up to 50e60% of capacity [7]. Maximal variation occurs in case of very high winds, i.e. when wind turbines are being stopped from rated to zero power production. The stopping of an entire wind farm usually does not happen instantaneously, i.e. it can last from several minutes to hours, depending on the wind farm size. Hourly wind power variations over the course of one day play an important role in the scheduling and dispatch of power plants and provision of secondary and tertiary reserves [8]. Thereby the impacts of wind power are dependent on forecast accuracy, as well as on generation mix and flexibility of the power system [9]. The long-term variability of wind power (for characteristic time scale from months to years) is more generally determined by seasonal meteorological patterns and by inter-annual variations of the wind. Since wind power plants are privileged in generation scheduling, their power is used whenever available. This means that high wind power penetration poses a great challenge for the transmission system operator, who is responsible for ensuring the reserves which are necessary to keep the system operating under a pre-stated range of conditions [10]. According to the EWEA report [9], when approximately 10% of total electricity consumption is produced by wind power, the increase in extra reserves calculated for various countries and regions is 2e4% of the installed wind power capacity, assuming state-of-the-art forecasting techniques. Hence, crucial parameters for large-scale integration of wind energy are the provisioning of necessary reserves and flexibility, geographical dispersion of wind power and the availability of forecasting techniques.
Direct generation of wind power is used since it appears to be more practical, leading to a reduced number of states and a lower order Markov chain at equal power data resolution [2]. 3.1. Markov chain basics A stochastic process is a Markov process if the probability of the present state depends only on the number of given states in the past, but not on any other past state. If such a process is discretized over time, it is called a Markov chain. The order of a Markov chain is equal to number of past states that are taken into account. Therefore, a first order Markov chain is a mathematical model of a discrete-time random variable, where the future state depends only of the present state of the process. The transition probability of moving from state i to state j is defined with [11]:
P Xtþ1 ¼ xj Xt ¼ xi ¼ pij
(2)
Since these probabilities are constant in time and divided by equal time steps, it is possible to formulate the transition probability matrix P for a first order Markov chain as follows [4]: Xtþ1 /
2
p12 p22 « pm2
p11 6 p21 P ¼ Xt Y6 4 « pm1
3 . p1m . p2m 7 7 : 1 « 5 . pmm
(3)
Matrix P is of size m m, where m is the number of defined states of variable X. The probabilities are obtained by computing all i,j transitions in the time series with:
nij pij ¼ Pm
j¼1
nij
;
0 pij 1
(4)
Thus, each row of matrix P represents the probability density function for the transition from a state i to a state j, where: m X
pij ¼ 1;
ci ¼ 1; .m
(5)
j¼1
For second-order Markov chain transition probabilities are defined as:
P Xtþ1 ¼ xk Xt1 ¼ xi ; Xt ¼ xj ¼ pijk ;
(6)
where pijk is the probability of the next state k if the current state is j and the previous state was i. In case of a second-order Markov chain, matrix P is symbolically shown as an m m2 sized matrix [4]:
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
2
p111 6 p121 6 6 « 6 6 p1m1 P ¼ 6 6 p211 6 6 p221 6 4 « pmm1 m X
p112 p122 « p1m2 p212 p222 « pmm2
/ / 1 / / / 1 /
3
p11m p12m 7 7 « 7 7 pi1m 7 7; p21m 7 7 p22m 7 7 « 5 pmmm
(7)
2605
using a speed/power curve. Due to the complexity of the terrain characteristics, the shading effect inside the wind farm and the different wind conditions at each wind turbine, a corresponding speed/power curve for the whole wind farm could be used, adjusted to include all these factors. If existing wind farm generation data is available, the best way is to take the total power output time series as the input data. 4.1. Resolution of the states discretization
pijk ¼ 1;
ci; j ¼ 1; .; m:
(8)
k¼1
3.2. MCMC simulation Markov chain Monte Carlo methods are a class of algorithms used for sampling from probability distributions based on Markov chain. Monte Carlo is a technique which involves using random numbers to obtain data series with a given probability. The first step in an MCMC simulation is to define the states of the associated Markov chain and the transition probability matrix. A higher number of defined states leads to better model accuracy, but significantly increases complexity. As opposed to modeling in the wind speed domain, when modeling in the wind power domain, there is significant state reduction since all the wind speeds in the zero and nominal power area correspond to these two states. A further simulation step requires construction of the cumulative probability transition matrix Pcum from the relation [2]:
pcum;ij ¼
j X
pil
(9)
To define an appropriate Markov chain and transition probability matrix P, wind power needs to be discretized into states. It is significant that higher values in matrix P are gathered around the main diagonal indicating high chronological persistence in the process. Namely, wind speed has a tendency to remain almost the same at subsequent moments in time. When defining the states, it is important to note that a higher resolution of discretization, i.e. a higher number of states, makes the model more complex, but more accurate. Since short-term variations are typically within 5e10% of the rated power, it follows that the width between the states should be significantly smaller than these variations in order to obtain good representation. Fig. 3 shows how the pdf of simulated 10 min wind power variations changes with different discretization of the states, where the width of the states is defined in relative terms to the installed power. In the performed simulations, a state width of 0.8% of installed capacity was used regardless of higher model complexity since it provided an excellent fit of the generated power variations pdf. On the other hand, if short-term power variations are not the main objective of the simulation, it is not necessary to use such a high discretization resolution.
l¼1
Each row i of Pcum corresponds to the discrete cumulative distribution function (cdf) for the next transition. The simulation starts with arbitrarily assumed initial wind power, i.e. initial state i. To sample the next state, random values between 0 and 1 are generated using a uniform random number generator. The value of the random number is compared with elements of the ith row of the Pcum matrix. Assuming this number falls between elements j 1 and j, then state j is chosen as the next state. The actual wind power inside state j can be obtained using the following relation [2]:
P ¼ Pj1 þ u Pj Pj1 ;
(10)
where Pj1 and Pj are the wind power boundaries of state j and u is the uniform random number. To sample the next transition, the same procedure is repeated for row j. In case of a second-order Markov chain, the simulation procedure is analogous, but there are two initial states assumed and the row index is defined by the current and preceding states.
4.2. The acf and wind power variations An autocorrelation function (acf) is generally used to determine the persistence structure of observed data. In the case of wind power data series, its acf can show how variations change on different time scales. Fig. 4 shows the acf of the original and two generated wind power data series. It can be seen that both generated data series (using 1st and 2nd order Markov chains) fail to retain the persistence of the original data on a longer time scale. Furthermore, the exponential fall of autocorrelation is faster for the first order Markov chain. Consequently, the acf of the generated data indicates that this model can adequately simulate only shortterm (10 min) wind power variations. For longer time-scales, the generated data loses persistence, creating wind power variations
4. Simulation settings In the experiments performed in this work, the MCMC simulation is used to generate wind power time series at two sites, with the main objective that the total power and its short-term (10 min) variations mimic the original data series as closely as possible. Assessment of the simulation results is based on: e the cdf of the total power, e the pdf of total power short-term variations. The dataset used was a time series of 10 min wind speed averages measured over one year (6 24 365 ¼ 52,560 data points) at two wind farm sites. These speeds are transformed to wind power
Fig. 3. Comparison of the wind power variations pdf for different states discretization.
2606
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
original data with respect to the statistical properties described in the previous section. To assess the quality of the generated data, the cumulative distribution function of total power and the probability density function of its 10 min variations are implemented. 5.1. Independent wind power series generation
Fig. 4. The autocorrelation functions of the original and generated wind power data series.
significantly higher than the original. Fig. 5 shows that the distribution of the generated data 10 min variations fits perfectly with original data, while the pdf of hourly variations shows significant misalignment. However, a representation of the acf with a closer fit of hourly variations could be obtained using a data of hourly averages, i.e. a 1 hour wind power averaging window. Since the original dataset available for this work was time series of 10 min averages, in the remainder of this work only 10 min variations are observed to make full use of the granularity of the samples, offering most representative available short term variations. However, longer term variations could be analyzed equivalently, but using a longer wind power averaging window, in order to obtain satisfactory reproduction of the acf [2]. 5. Preliminary simulation analysis To successfully analyze the effects of varying wind power production on a power system, it is necessary to consider power output from multiple wind farms in a relevant geographical area. This additionally complicates the simulation model since it requires modeling of the stochastic dependence between each wind farm power production. In this paper, the annual measurements of wind data from two wind farm sites are used. The distance between the sites is about 90 km and the linear correlation coefficient between their power output is r ¼ 0.65. The MCMC simulation is used to synthetically generate their power output, with the main objective to mimic the
Fig. 5. Comparison of the pdf for 10 min and hourly variations of the generated and original wind power.
The first simulation is performed independently for each of the two wind farms. This means that a first order Markov chain is defined separately for each site, with independent transition probability matrices and uniform random number generators. This simulation method gives two data series whose linear correlation coefficient is approximately zero. Alternatively, if the same random generator is used (u2 ¼ u1), the linear correlation coefficient between the two generated wind power time series rises to around rz0:9. Additionally, if the second random number generator is set to u2 ¼ 1 u1, it results with a negative linear correlation coefficient of rz 0:4. These three simulations are compared with the original data with respect to the cdf of total power distribution (Fig. 6) and the pdf of its 10 min variations (Fig. 7). In case of rz0, the independency between the two series gives a lower probability of higher total power rates and a higher probability of lower power rates, i.e. deviation of the cdf shape, in comparison with the original one. A negative linear correlation value of rz 0:4 causes even more frequent lower power rates, resulting in total power variations smaller then the original (Fig. 7). With a correlation coefficient of about rz0:9, the simulated data series is more dependent than the original one. Consequently, the cdf is not approximated well in this case either (Fig. 6). It is important to note that considering more wind farm sites in the simulation, mismatch of the simulated total power variations pdf (Fig 7) would be even more significant, since more dependency relations are ignored using this simulation method. These simulation results indicate that assuming independence when modeling power production of multiple wind farms can lead to certain fallacies. Hence, the linear correlation of the generated wind power needs to be as close as possible to the original (r ¼ 0.65), in order to obtain the correct cdf of the total power production. Furthermore, considering linear correlation as the only measure of stochastic dependence may lead to the false conclusions. Linear correlation can be used to show the dependency between normally distributed variables, but in cases of non-normal distributions, it cannot capture the non-linear, monotonic
Fig. 6. The cdf for the original and generated total wind power with different linear correlations.
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
2607
Fig. 8. The empirical copula of the two observed sites. Fig. 7. The pdf for 10 min variations of the original and generated total wind power variations.
relationships between variables. In order to capture such a relationship, the one-dimensional marginal distribution should be transformed into ranks, by applying the cumulative distribution function. Generally, the modeling procedure is decomposed into two basic steps, i.e. modeling of the marginal distributions (cdf and pdf) and modeling of the stochastic dependence structure. Based on this idea, specific functions, called copula functions, can be used to model the dependence between two variables more accurately [12]. 5.2. Wind power series generation using copulas A copula is a multivariate joint distribution defined on an ndimensional unit cube [0,1]n such that every marginal distribution is uniform over the interval [0,1]. Specifically, for the bivariate case, Sklar’s theorem says: let (X,Y) be a bivariate random vector with a joint cdf Fxy(X,Y) and marginal cdfs Fx(x) and Fy(y). Then there exists a copula C such that [13]:
FXY ðx; yÞ ¼ CðFX ðxÞ; FY ðyÞÞ
(11)
There are many families of copulas which differ with respect to the aspect of dependence they represent. In this work, the following empirical copula is used: Let (Xi,Yi), i ¼ 1,.,n be a sample from (X,Y). Then the empirical distribution function is:
Fn ðx; yÞ ¼
n 1X IfXi x; Yi yg; n i¼1
used to describe the dependency between the two generated time series. Simulation for the first site is performed using the same MCMC simulation described previously, while the second site wind power series is obtained from the empirical copula using a Monte Carlo simulation. The linear correlation coefficient between the wind power series generated by this simulation was equal to the original one (r ¼ 0.65), which gave a cdf of total wind power which almost perfectly describes the original (Fig. 9). However, Fig. 10 shows that the pdf of these variations is dramatically off. The reason for this is the fact that the copula describes the dependence between the two time series, but chronological persistence inside the second time series is ignored. Thus, the power variations at the second site appear to be abnormally high. The simulation results obtained for independent first-order MCMC simulation indicate that the total power distribution can not be correctly generated without modeling the stochastic dependence between the productions at each wind farm. In this paper, using copulas to describe their stochastic dependence is proposed. This approach has been shown to obtain good solution for the total power distribution and energy capacity, but fails to adequately model the associated power variations. 6. The proposed modified second-order MCMC approach and simulation results Herein, a modified MCMC simulation approach is presented aimed to generate multiple time series while preserving their
(12)
and marginal cdfs are:
Fnx ðxÞ ¼ Fn ðx; þNÞ; Fny ðyÞ ¼ Fn ð þ N; yÞ;
(13)
where: N < x; y < þN. The empirical copula Cn is defined by:
1 1 Cn ðu; vÞ ¼ Fn Fnx ðuÞ; Fny ðvÞ :
(14)
For the wind power time series at the two considered sites, the empirical copula is calculated and graphically presented using a scatter diagram (Fig. 8). It is evident that a higher concentration occurs at zero and nominal output power. Furthermore, there is a low probability of high wind speed at wind farm site 1 at the same time as low wind speed at site 2. Hence, the empirical copula is
Fig. 9. The cdf of the original and copula-generated total wind power.
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
2608 original series
copula-generated
0,7
0,6
Probability
0,5
0,4
0,3
0,2
0,1
0,0 -10%
-5%
0%
5%
10%
Total wind power variations [%of capacity]
Fig. 10. The pdf for 10 min variations of the original and copula-generated total wind power.
individual persistence, as well as their mutual dependency. All the considered time series are modeled as one Markov process, but with separate transition probability matrices for each series. It is important to note that transition probabilities are defined in such a way that the future state of each time series depends on the previous state and recent states in all other time series. Therefore, the order of the Markov chain and number of transition probability matrices needs to be at least equal to the number of the considered time series. The MCMC simulation is preformed in the same way as described in Section 3, but using a different cumulative probability matrix for each time series. For better understanding of the simulation algorithm, a worked example using the two considered real wind farm sites is described and tested. These two wind power data series are considered as one second-order Markov chain, while transitions between the states are described by two matrices. The transition probabilities of the first matrix Px define the state appearance in the first time series:
o n P Xtþ1 ¼ xk Xt ¼ xi ; Yt ¼ yj ¼ pxijk ;
(15)
while the second one Py in the second time series:
o n P Ytþ1 ¼ yl Xtþ1 ¼ xk ; Yt ¼ yj ¼ pyjkl :
Fig. 12. The cdf for the original and proposed 2nd order MCMC generated total wind power.
Furthermore, a future state in the second series depends on the present state in first time series and previous state in second time series. Using the MCMC simulation, wind powers are generated interchangeably for both series, using matrix Px for the first site and Py for the second one. The transitions and the simulation flow are illustrated in Fig. 11. This simulation method provides significant improvement in the resulting aggregate wind power data series. In Fig. 12, it is shown that the cdf of the simulated total wind power perfectly matches the original one, with a correlation coefficient fully corresponding to the original (r ¼ 0.65). In comparison to the MCMC simulation with copula implementation, this method also obtains superior results for the pdf of 10 min variations shown in Fig. 13. Hence, unlike other simulation methods analyzed here, using the proposed MCMC approach it is possible to generate wind data almost fully corresponding to the original one, allowing satisfactory representation of short term wind power variations inside both time series, along with their aggregate power. This simulation method can generally be analogously applied to more wind farm sites, where the order of the Markov chain needs to be at least equal to number of considered sites. Transition probabilities for n matrices can be written as:
(16)
Hence, a future state in the first time series depends on the previous states in both the first and second time series.
Probability
0,7
0,6
proposed 2nd order MCMC original series
0,5
independent MCMC (ρ=0) copula generated
0,4
0,3
0,2
0,1
0,0 -10%
-5%
0%
5%
10%
Total wind power variations [% of capacity]
Fig. 11. Illustration of the proposed second-order MCMC simulation.
Fig. 13. The pdf for 10 min variations of the original and proposed 2nd order MCMC generated total wind power.
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
1st site x1-i
t
px1-ij...kl
2nd site x2-j
7. Conclusion
nth site
...
xn-k
px2-j...klm x1-l
t+1
...
x1-m
.. .
xn-o
pxn-klm...o
.. .
.. .
Fig. 14. Illustration of the proposed nth-order MCMC simulation.
P X1
tþ1
¼ x1l X1
¼ x1i ; X2
t
t
¼ x2j ; .; Xn
t
¼ xnk
¼ px1ij.kl (17) P X2
tþ1
¼ x2m X2
t
¼ x2j ; .; Xn
t
¼ xnk ; X1
tþ1
¼ x1l
¼ px2j.klm (18)
.
PfXn
tþ1
¼ xno jXn
t
¼ xnk ; X1
¼ x2m ; .g ¼ pxn
tþ1
klm.o
¼ x1m ; X2
2609
Analysis of the impact of wind power variations on power system operation requires adequate modeling of total power output from geographically dispersed wind farms. In this paper, different simulation methods are investigated on annual 10 min wind data series for two real wind sites. The obtained results show that adequate stochastic dependence modeling is required in order to obtain correct distributions of aggregate power output and total power variations. Additionally, for fitting the pdf of short-term variations, it is necessary to preserve at least short-term persistence inside generated power time series. A modified second-order Markov chain Monte Carlo simulation method is proposed in order to obtain synthetic time series of aggregate wind power which fits the original data very well with respect to the basic statistical parameters. This method can show certain limitations in considerable memory usage based on the number of considered wind sites, but does not play a major role in the preprocessing and main algorithm complexity. For future work it is planned to additionally investigate the possibility of using this method for aggregate wind farms generation modeling. Furthermore, it will be assessed whether the proposed method can also be used to model the dependency between wind power generation and load.
tþ1
(19)
The mutual dependence between wind data series at multiple sites could be simply examined using a correlation matrix containing all combinations of cross-correlations. Illustration of the simulation algorithm for n sites is shown in Fig. 14. The main limitation of the presented simulation method lies in potentially high memory space requirements. Namely, since the order of the Markov chain is equal to the number of time series, the size of the transition probability matrix exponentially grows with number of wind farms sites. Since high resolution of the states is required (as described in Subsection 4.1), the dimensions of the transition matrices can be very large, even for a relatively low number of sites. However, the complexity of the simulation approach is equal to O(m n), where m is the number of states, and n is the number of values in the data series. This complexity also refers to the preprocessing transition matrices generation step. Consequently, the number of considered sites does not play a major role in the execution time of the preprocessing and the simulation itself, although it does significantly increase the required memory usage. Furthermore, the size of the transition probability matrix, i.e. the memory usage, can significantly be reduced if short-term wind power variations are not considered the main simulation objective. It follows that the granularity of the states defined for the considered number of sites, creates a trade-off between simulation complexity and accuracy.
References [1] World Wind Energy conference and exhibition. World wind energy report 2008, www.wwindea.org/home/images/stories/worldwindenergyreport2008 _s.pdf; June 2009. [2] Papaefthymiou G, Klöckl B. MCMC for wind power simulation. IEEE Transactions on Energy Conversion 2008;23(1). [3] Ettoumi FY, Sauvageot H, Adane A-E-H. Statistical bivariate modeling of wind using first-order Markov chain and Weibull distribution. Renewable Energy 2003;28:1787e802. [4] Shamshad A, Bawadi MA, Hussin W, Majid T, Sanusi S. First and second order Markov chain models for synthetic generation of wind speed time series. Energy 2005;30:693e708. [5] Sahin AD, Sen Z. First-order Markov chain approach to wind speed modeling. Journal of Wind Engineering and Industrial Aerodynamics 2001;89:263e9. [6] Van der Hoven I. Power spectrum of horizontal wind speed in the frequency range from 0.0007 to 900 cycles per hour. Journal of Meteorology 1956;14. [7] Nørgård P, Giebel G, Holttinen H, Söder L, Petterteig A. Fluctuations and predictability of wind and hydropower, Risø-R-1443, 2004. [8] Goic R, Jakus D, Penovic I, Vasilj J, Kovacevic A, Lovric M, Klaric M. Integration of wind energy into the Croatian power system, Croatian national committee of CIGRE, Symposium on power system management, Cavtat 2008. [9] Holttinen H. Technical research Centre of Finland VTT, Impact of hourly wind Power variations on the system operation in Nordic countries, P.O. Box 1601, FIN-02044 VTT, 2005. [10] European wind energy association (EWEA). Large scale integration of wind energy in the European power supply: analysis, issues and recommendations, EWEA report, http://www.ewea.org/fileadmin/ewea_documents/documents/ publications/grid/051215_Grid_report.pdf; 2005. [11] Meyn SP, Tweedie RL. Markov chains and stochastic stability. first ed. London: Springer-Verlag; 1993. [12] Papaefthymiou G. Integration of stochastic generation in power systems, Ph.D. dissertation, Delft university of technology, ISBN: 978-90-8570-186-6; 2007. [13] Fermanian J, Radulovic D, Wegkamp M. Weak convergence of empirical copula processes. Bernoulli 2004;10(5):847e60.
Contents lists available at ScienceDirect
Renewable Energy journal homepage: www.elsevier.com/locate/renene
Simulation of aggregate wind farm short-term production variations R. Goi c, J. Krstulovi c*, D. Jakus Faculty of Electrical Engineering and Naval Architecture, University of Split, Rudera Boskovica bb, 21000 Split, Croatia
a r t i c l e i n f o
a b s t r a c t
Article history: Received 25 November 2009 Accepted 7 April 2010 Available online 1 May 2010
The variability of wind power production poses the greatest challenge in the integration of large-scale wind power in power systems. Furthermore, larger-scale penetration implies a wider geographical spreading of installed wind power, resulting in reduced variability and the smoothing effect of total power generation. Therefore, analysis of the impact of wind power variations on power system operation requires adequate modeling of aggregate power output from geographically dispersed wind farms. This paper analyzes different aspects of Markov chain Monte Carlo simulation methods for the synthetic generation of dependent wind power time series. However, testing indicates that these approaches do not adequately model the stochastic dependence between wind power time series in conjunction with individual persistence which is necessary to obtain realistic distributions of aggregate power output and total power variations. Consequently, a novel approach based on a modified second-order Markov chain Monte Carlo simulation is proposed. Simulation results show that this method obtains synthetic time series of aggregate wind power which very closely fit the original data, with respect to both the cumulative density function of total output power and the probability density function of power variations. Ó 2010 Elsevier Ltd. All rights reserved.
Keywords: Aggregate wind power generation Wind power variations Monte Carlo simulation Markov chain Stochastic dependency
1. Introduction With an installed nominal capacity of 121.2 GW, wind power makes about 1.5% of world wide electricity generation (up to the end of 2008). This number is growing rapidly, having doubled in the three years between 2005 and 2008. Several countries have achieved relatively high levels of wind penetration, such as 19% of electricity production in Denmark, 11% in Spain and Portugal, and 7% in Germany in 2008 [1]. With this trend of growing wind power capacities, the integration of large-scale wind power in power systems has become more difficult. The stochastic nature of wind creates fluctuating and intermitting wind power. Since instantaneous electrical generation and consumption must remain balanced to maintain grid stability, such fluctuations create substantial obstacles in incorporating large amounts of wind power into a grid. Variations of wind power production occur on all timescales, from seconds to years, where even short-term variations are, to some extent, unpredictable. Observing the power variations at a single wind turbine or a single wind farm may lead to false conclusions. Namely, the geographical spreading of wind power reduces variability, increases predictability and reduces occasions with zero or peak output.
* Corresponding author. Tel.: þ385 21 305 802. E-mail address: [email protected] (J. Krstulovi c). 0960-1481/$ e see front matter Ó 2010 Elsevier Ltd. All rights reserved. doi:10.1016/j.renene.2010.04.005
Consequently, as wind farms are spread over wider geographical areas working under different wind conditions, the resulting total wind power generation becomes smoother. Designing a wind farm requires installing a meteorological station at a potential location and recording all the relevant wind data series. These time series are usually not long enough to create a representative model for potential wind farm generation and may lead to biased conclusions regarding energy capacity estimation, cost effectiveness, system design issues, power variations, etc. Consequently, it is necessary to use stochastic simulation techniques to synthetically generate wind speed (wind power) data while preserving its statistical parameters. Since wind speed is conditioned with lots of climatologic, atmospheric and environmental parameters, such phenomena are impossible to describe sufficiently using a deterministic approach. Furthermore, when studying the effects of varying wind power production on the power system, it is important to consider larger areas, i.e. aggregate power output of multiple wind farms. This additionally complicates the problem since simulating wind power for two or more sites also requires modeling their stochastic dependence in order to obtain realistic total power generation data and its variations. In related literature, various probabilistic and stochastic methods are used for synthetic generation of time series. The most relevant methods include: Monte Carlo simulation of independent random numbers distributed according to a Weibull probability distribution function, autoregressive (AR) in combination with
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
moving average (ARMA) models, Markov chain Monte Carlo (MCMC) methods, and a Wavelet e based approach. In the case when a certain persistency is required in a generated time series, the most commonly used methods are ARMA and MCMC. The basic advantage of the AR and ARMA models is controlled synthetization of a desired autocorrelation function (acf) while defining a relatively small number of parameters. On the other hand, the probability density function (pdf) of a time series generated with ARMA, generally does not match that of the measured data, leading to a wrong estimation of the energy yield [2]. Appropriate use of the MCMC method for wind speed (power) time series simulation leads to satisfactory pdf and acf matching between the generated and measured data. In this method, wind activity is modeled as a Markov chain, i.e. a discrete stochastic process with defined probabilities for transitions between states. In recent works [2e5], different aspects of using the MCMC method for the simulation of wind speed data are analyzed. G. Papaefthymiou and B. Klöckl [2] show that obtaining a stochastic model directly in the wind power domain is superior to indirect modeling of wind power using wind speeds. This leads to a reduced number of states requiring a lower order Markov chain at equal power data resolution. In this paper, emphasis is given to the simulation of aggregate wind farm power output with correct representation of their total power and its short term (10 min) variations. Wind power time series are independently generated for two potential wind farm sites using the MCMC simulation method. Analysis of the aggregate power output distribution and short-term variations indicates that certain dependency between the generated time series is required to obtain satisfactory results. Furthermore, to describe their dependency, an empirical copula is implemented in the MCMC simulation, resulting in good representation of aggregate wind power, but lacking the variations. To overcome these drawbacks, an algorithm is proposed using a modified MCMC simulation for synthetic generation of aggregate wind farm power output. Simulation results show that this method provides realistic representation of both aggregate power and its short term variations distribution. All the simulation methods were tested using 10 min average wind speed data over one year from two potential wind farm sites in Croatia. 2. Wind power variations The stochastic nature of wind causes significantly variable and intermittent wind turbine power output. This variability is present on all timescales: seconds, minutes, hours, days, months and years. Successfully defining and forecasting such variations makes crucial step for wind power integration and achieving optimal operation of a power system. The variability of the horizontal winds, that wind turbines interact with, shows interesting general characteristics approved with their kinetic energy of the speed fluctuation distribution in the frequency domain, known as the van der Hoven spectrum. Van der Hoven [6] shows that independently of the site, the typical spectrum exhibits two peaks approximately at 0.01 cycles/ hour (4 days period), and 50 cycles/h (1 min period), which are separated by an energy gap between periods of 10 min and 2 h. In a 4 days period, the peak corresponds to wind speed fluctuations due to the migratory pressure systems of a synoptic weather-map scale. The other excerpted peak is in the micrometeorological range as a result of the turbulence associated to local winds. The reason for the spectral gap is in the lack of a physical process which could support wind speed fluctuations in that frequency range. When comparing wind fluctuations of 10 min and hourly
2603
averages, both are placed in the spectral gap, whereas energy is higher in the 10 min period, closer to the high frequency peak, while in the hourly period the spectral gap mostly reaches its minimum. In this paper, only short-term power variations are analyzed, which include timescales from minutes to hours. Second-tosecond variations of individual wind turbine power output are insignificant due to the inertia of large rotating blades and the varying speed of the rotor. If the total power generated by a single wind farm is considered as a whole, variations of several seconds become absorbed since the same gusts do not occur simultaneously at all turbines. Furthermore, the total wind power variations of multiple wind farms are not significant, even for a few minutes period, depending on their geographical spreading. Namely, in larger geographical areas, weather fronts can not pass through all wind farms simultaneously. Thus, longer time scales are subject to the so-called smoothing effect, as illustrated in Fig. 1. Consequently, the maximum amplitudes of wind power variations experienced in the power system are reduced. Fig. 2 shows that higher variations of 10 min wind power occur less frequently when the total power of three wind farms is considered, in comparison to that of a single wind farm. This effect would be much more significant for a higher number of wind farms over an even larger geographic area. For better understanding of the geographical dispersion effect on short term wind power variations, it is important to observe the linear correlation coefficient between two wind power time series. Namely, dependence is most often measured with a linear correlation coefficient. This coefficient generally indicates the strength and direction of a linear relationship between two random variables and is calculated by:
P
ðx xÞðy yÞ
ffi rx;y ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi P P ðx xÞ2
ðy yÞ2
(1)
The linear correlation coefficient can range from 1 to 1, where the two extreme cases and the middle value are the most interesting to observe. For a maximal value of the correlation coefficient r ¼ 1, increases and decreases in power production from individual outputs occur simultaneously, resulting in the largest possible variations of total power production. For the case where r ¼ 0, the power outputs fluctuate independently, causing smaller variations. For the third case, when the correlation coefficient is r ¼ 1, the outputs fluctuate in opposite directions, resulting in minimal total
Fig. 1. An example of the smoothing effect by geographic dispersion for three wind farm locations in Croatia.
2604
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
3. Markov chain Monte Carlo simulation Available data series recorded at real wind farm locations are generally not long enough to be representative and sufficient for a detailed analysis. Thus, various simulation techniques are used for synthetic generation of new data. In this paper, the Markov chain Monte Carlo simulation method is applied. Using this approach, wind speed data is modeled as a Markov chain, while the Monte Carlo simulation is used for generating new data. The main objective of this simulation is to obtain realistic wind power output time series and can be obtained using one of two approaches: e first building a model in the wind speed domain and then transforming generated wind speeds to wind power or, e first transforming available wind speed data into wind power and then directly generating wind power time series. Fig. 2. The probability density function of 10 min wind power variations e the effect of geographical dispersion.
power variations. Therefore, smaller variations are obtained for lower correlation coefficients. In reality, this coefficient ranges between zero and one, depending on the geographical distance between the observed sites with inverse exponential dependency. Furthermore, considering multiple wind farm sites, correlation between sites could be represented using a correlation square matrix whose size is equal to number of sites. Each row of the correlation matrix represents the correlation between one site and all the other sites. Thus, the diagonal is composed of all ones while the remaining elements contain all cross-correlations between sites. Generally, fast power variations (within a minute time interval) are usually very low, with no impact on power system operation. On the other hand, in-hour variations (10e60 min) can be much more significant, with extreme ramp rates reaching up to 50e60% of capacity [7]. Maximal variation occurs in case of very high winds, i.e. when wind turbines are being stopped from rated to zero power production. The stopping of an entire wind farm usually does not happen instantaneously, i.e. it can last from several minutes to hours, depending on the wind farm size. Hourly wind power variations over the course of one day play an important role in the scheduling and dispatch of power plants and provision of secondary and tertiary reserves [8]. Thereby the impacts of wind power are dependent on forecast accuracy, as well as on generation mix and flexibility of the power system [9]. The long-term variability of wind power (for characteristic time scale from months to years) is more generally determined by seasonal meteorological patterns and by inter-annual variations of the wind. Since wind power plants are privileged in generation scheduling, their power is used whenever available. This means that high wind power penetration poses a great challenge for the transmission system operator, who is responsible for ensuring the reserves which are necessary to keep the system operating under a pre-stated range of conditions [10]. According to the EWEA report [9], when approximately 10% of total electricity consumption is produced by wind power, the increase in extra reserves calculated for various countries and regions is 2e4% of the installed wind power capacity, assuming state-of-the-art forecasting techniques. Hence, crucial parameters for large-scale integration of wind energy are the provisioning of necessary reserves and flexibility, geographical dispersion of wind power and the availability of forecasting techniques.
Direct generation of wind power is used since it appears to be more practical, leading to a reduced number of states and a lower order Markov chain at equal power data resolution [2]. 3.1. Markov chain basics A stochastic process is a Markov process if the probability of the present state depends only on the number of given states in the past, but not on any other past state. If such a process is discretized over time, it is called a Markov chain. The order of a Markov chain is equal to number of past states that are taken into account. Therefore, a first order Markov chain is a mathematical model of a discrete-time random variable, where the future state depends only of the present state of the process. The transition probability of moving from state i to state j is defined with [11]:
P Xtþ1 ¼ xj Xt ¼ xi ¼ pij
(2)
Since these probabilities are constant in time and divided by equal time steps, it is possible to formulate the transition probability matrix P for a first order Markov chain as follows [4]: Xtþ1 /
2
p12 p22 « pm2
p11 6 p21 P ¼ Xt Y6 4 « pm1
3 . p1m . p2m 7 7 : 1 « 5 . pmm
(3)
Matrix P is of size m m, where m is the number of defined states of variable X. The probabilities are obtained by computing all i,j transitions in the time series with:
nij pij ¼ Pm
j¼1
nij
;
0 pij 1
(4)
Thus, each row of matrix P represents the probability density function for the transition from a state i to a state j, where: m X
pij ¼ 1;
ci ¼ 1; .m
(5)
j¼1
For second-order Markov chain transition probabilities are defined as:
P Xtþ1 ¼ xk Xt1 ¼ xi ; Xt ¼ xj ¼ pijk ;
(6)
where pijk is the probability of the next state k if the current state is j and the previous state was i. In case of a second-order Markov chain, matrix P is symbolically shown as an m m2 sized matrix [4]:
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
2
p111 6 p121 6 6 « 6 6 p1m1 P ¼ 6 6 p211 6 6 p221 6 4 « pmm1 m X
p112 p122 « p1m2 p212 p222 « pmm2
/ / 1 / / / 1 /
3
p11m p12m 7 7 « 7 7 pi1m 7 7; p21m 7 7 p22m 7 7 « 5 pmmm
(7)
2605
using a speed/power curve. Due to the complexity of the terrain characteristics, the shading effect inside the wind farm and the different wind conditions at each wind turbine, a corresponding speed/power curve for the whole wind farm could be used, adjusted to include all these factors. If existing wind farm generation data is available, the best way is to take the total power output time series as the input data. 4.1. Resolution of the states discretization
pijk ¼ 1;
ci; j ¼ 1; .; m:
(8)
k¼1
3.2. MCMC simulation Markov chain Monte Carlo methods are a class of algorithms used for sampling from probability distributions based on Markov chain. Monte Carlo is a technique which involves using random numbers to obtain data series with a given probability. The first step in an MCMC simulation is to define the states of the associated Markov chain and the transition probability matrix. A higher number of defined states leads to better model accuracy, but significantly increases complexity. As opposed to modeling in the wind speed domain, when modeling in the wind power domain, there is significant state reduction since all the wind speeds in the zero and nominal power area correspond to these two states. A further simulation step requires construction of the cumulative probability transition matrix Pcum from the relation [2]:
pcum;ij ¼
j X
pil
(9)
To define an appropriate Markov chain and transition probability matrix P, wind power needs to be discretized into states. It is significant that higher values in matrix P are gathered around the main diagonal indicating high chronological persistence in the process. Namely, wind speed has a tendency to remain almost the same at subsequent moments in time. When defining the states, it is important to note that a higher resolution of discretization, i.e. a higher number of states, makes the model more complex, but more accurate. Since short-term variations are typically within 5e10% of the rated power, it follows that the width between the states should be significantly smaller than these variations in order to obtain good representation. Fig. 3 shows how the pdf of simulated 10 min wind power variations changes with different discretization of the states, where the width of the states is defined in relative terms to the installed power. In the performed simulations, a state width of 0.8% of installed capacity was used regardless of higher model complexity since it provided an excellent fit of the generated power variations pdf. On the other hand, if short-term power variations are not the main objective of the simulation, it is not necessary to use such a high discretization resolution.
l¼1
Each row i of Pcum corresponds to the discrete cumulative distribution function (cdf) for the next transition. The simulation starts with arbitrarily assumed initial wind power, i.e. initial state i. To sample the next state, random values between 0 and 1 are generated using a uniform random number generator. The value of the random number is compared with elements of the ith row of the Pcum matrix. Assuming this number falls between elements j 1 and j, then state j is chosen as the next state. The actual wind power inside state j can be obtained using the following relation [2]:
P ¼ Pj1 þ u Pj Pj1 ;
(10)
where Pj1 and Pj are the wind power boundaries of state j and u is the uniform random number. To sample the next transition, the same procedure is repeated for row j. In case of a second-order Markov chain, the simulation procedure is analogous, but there are two initial states assumed and the row index is defined by the current and preceding states.
4.2. The acf and wind power variations An autocorrelation function (acf) is generally used to determine the persistence structure of observed data. In the case of wind power data series, its acf can show how variations change on different time scales. Fig. 4 shows the acf of the original and two generated wind power data series. It can be seen that both generated data series (using 1st and 2nd order Markov chains) fail to retain the persistence of the original data on a longer time scale. Furthermore, the exponential fall of autocorrelation is faster for the first order Markov chain. Consequently, the acf of the generated data indicates that this model can adequately simulate only shortterm (10 min) wind power variations. For longer time-scales, the generated data loses persistence, creating wind power variations
4. Simulation settings In the experiments performed in this work, the MCMC simulation is used to generate wind power time series at two sites, with the main objective that the total power and its short-term (10 min) variations mimic the original data series as closely as possible. Assessment of the simulation results is based on: e the cdf of the total power, e the pdf of total power short-term variations. The dataset used was a time series of 10 min wind speed averages measured over one year (6 24 365 ¼ 52,560 data points) at two wind farm sites. These speeds are transformed to wind power
Fig. 3. Comparison of the wind power variations pdf for different states discretization.
2606
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
original data with respect to the statistical properties described in the previous section. To assess the quality of the generated data, the cumulative distribution function of total power and the probability density function of its 10 min variations are implemented. 5.1. Independent wind power series generation
Fig. 4. The autocorrelation functions of the original and generated wind power data series.
significantly higher than the original. Fig. 5 shows that the distribution of the generated data 10 min variations fits perfectly with original data, while the pdf of hourly variations shows significant misalignment. However, a representation of the acf with a closer fit of hourly variations could be obtained using a data of hourly averages, i.e. a 1 hour wind power averaging window. Since the original dataset available for this work was time series of 10 min averages, in the remainder of this work only 10 min variations are observed to make full use of the granularity of the samples, offering most representative available short term variations. However, longer term variations could be analyzed equivalently, but using a longer wind power averaging window, in order to obtain satisfactory reproduction of the acf [2]. 5. Preliminary simulation analysis To successfully analyze the effects of varying wind power production on a power system, it is necessary to consider power output from multiple wind farms in a relevant geographical area. This additionally complicates the simulation model since it requires modeling of the stochastic dependence between each wind farm power production. In this paper, the annual measurements of wind data from two wind farm sites are used. The distance between the sites is about 90 km and the linear correlation coefficient between their power output is r ¼ 0.65. The MCMC simulation is used to synthetically generate their power output, with the main objective to mimic the
Fig. 5. Comparison of the pdf for 10 min and hourly variations of the generated and original wind power.
The first simulation is performed independently for each of the two wind farms. This means that a first order Markov chain is defined separately for each site, with independent transition probability matrices and uniform random number generators. This simulation method gives two data series whose linear correlation coefficient is approximately zero. Alternatively, if the same random generator is used (u2 ¼ u1), the linear correlation coefficient between the two generated wind power time series rises to around rz0:9. Additionally, if the second random number generator is set to u2 ¼ 1 u1, it results with a negative linear correlation coefficient of rz 0:4. These three simulations are compared with the original data with respect to the cdf of total power distribution (Fig. 6) and the pdf of its 10 min variations (Fig. 7). In case of rz0, the independency between the two series gives a lower probability of higher total power rates and a higher probability of lower power rates, i.e. deviation of the cdf shape, in comparison with the original one. A negative linear correlation value of rz 0:4 causes even more frequent lower power rates, resulting in total power variations smaller then the original (Fig. 7). With a correlation coefficient of about rz0:9, the simulated data series is more dependent than the original one. Consequently, the cdf is not approximated well in this case either (Fig. 6). It is important to note that considering more wind farm sites in the simulation, mismatch of the simulated total power variations pdf (Fig 7) would be even more significant, since more dependency relations are ignored using this simulation method. These simulation results indicate that assuming independence when modeling power production of multiple wind farms can lead to certain fallacies. Hence, the linear correlation of the generated wind power needs to be as close as possible to the original (r ¼ 0.65), in order to obtain the correct cdf of the total power production. Furthermore, considering linear correlation as the only measure of stochastic dependence may lead to the false conclusions. Linear correlation can be used to show the dependency between normally distributed variables, but in cases of non-normal distributions, it cannot capture the non-linear, monotonic
Fig. 6. The cdf for the original and generated total wind power with different linear correlations.
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
2607
Fig. 8. The empirical copula of the two observed sites. Fig. 7. The pdf for 10 min variations of the original and generated total wind power variations.
relationships between variables. In order to capture such a relationship, the one-dimensional marginal distribution should be transformed into ranks, by applying the cumulative distribution function. Generally, the modeling procedure is decomposed into two basic steps, i.e. modeling of the marginal distributions (cdf and pdf) and modeling of the stochastic dependence structure. Based on this idea, specific functions, called copula functions, can be used to model the dependence between two variables more accurately [12]. 5.2. Wind power series generation using copulas A copula is a multivariate joint distribution defined on an ndimensional unit cube [0,1]n such that every marginal distribution is uniform over the interval [0,1]. Specifically, for the bivariate case, Sklar’s theorem says: let (X,Y) be a bivariate random vector with a joint cdf Fxy(X,Y) and marginal cdfs Fx(x) and Fy(y). Then there exists a copula C such that [13]:
FXY ðx; yÞ ¼ CðFX ðxÞ; FY ðyÞÞ
(11)
There are many families of copulas which differ with respect to the aspect of dependence they represent. In this work, the following empirical copula is used: Let (Xi,Yi), i ¼ 1,.,n be a sample from (X,Y). Then the empirical distribution function is:
Fn ðx; yÞ ¼
n 1X IfXi x; Yi yg; n i¼1
used to describe the dependency between the two generated time series. Simulation for the first site is performed using the same MCMC simulation described previously, while the second site wind power series is obtained from the empirical copula using a Monte Carlo simulation. The linear correlation coefficient between the wind power series generated by this simulation was equal to the original one (r ¼ 0.65), which gave a cdf of total wind power which almost perfectly describes the original (Fig. 9). However, Fig. 10 shows that the pdf of these variations is dramatically off. The reason for this is the fact that the copula describes the dependence between the two time series, but chronological persistence inside the second time series is ignored. Thus, the power variations at the second site appear to be abnormally high. The simulation results obtained for independent first-order MCMC simulation indicate that the total power distribution can not be correctly generated without modeling the stochastic dependence between the productions at each wind farm. In this paper, using copulas to describe their stochastic dependence is proposed. This approach has been shown to obtain good solution for the total power distribution and energy capacity, but fails to adequately model the associated power variations. 6. The proposed modified second-order MCMC approach and simulation results Herein, a modified MCMC simulation approach is presented aimed to generate multiple time series while preserving their
(12)
and marginal cdfs are:
Fnx ðxÞ ¼ Fn ðx; þNÞ; Fny ðyÞ ¼ Fn ð þ N; yÞ;
(13)
where: N < x; y < þN. The empirical copula Cn is defined by:
1 1 Cn ðu; vÞ ¼ Fn Fnx ðuÞ; Fny ðvÞ :
(14)
For the wind power time series at the two considered sites, the empirical copula is calculated and graphically presented using a scatter diagram (Fig. 8). It is evident that a higher concentration occurs at zero and nominal output power. Furthermore, there is a low probability of high wind speed at wind farm site 1 at the same time as low wind speed at site 2. Hence, the empirical copula is
Fig. 9. The cdf of the original and copula-generated total wind power.
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
2608 original series
copula-generated
0,7
0,6
Probability
0,5
0,4
0,3
0,2
0,1
0,0 -10%
-5%
0%
5%
10%
Total wind power variations [%of capacity]
Fig. 10. The pdf for 10 min variations of the original and copula-generated total wind power.
individual persistence, as well as their mutual dependency. All the considered time series are modeled as one Markov process, but with separate transition probability matrices for each series. It is important to note that transition probabilities are defined in such a way that the future state of each time series depends on the previous state and recent states in all other time series. Therefore, the order of the Markov chain and number of transition probability matrices needs to be at least equal to the number of the considered time series. The MCMC simulation is preformed in the same way as described in Section 3, but using a different cumulative probability matrix for each time series. For better understanding of the simulation algorithm, a worked example using the two considered real wind farm sites is described and tested. These two wind power data series are considered as one second-order Markov chain, while transitions between the states are described by two matrices. The transition probabilities of the first matrix Px define the state appearance in the first time series:
o n P Xtþ1 ¼ xk Xt ¼ xi ; Yt ¼ yj ¼ pxijk ;
(15)
while the second one Py in the second time series:
o n P Ytþ1 ¼ yl Xtþ1 ¼ xk ; Yt ¼ yj ¼ pyjkl :
Fig. 12. The cdf for the original and proposed 2nd order MCMC generated total wind power.
Furthermore, a future state in the second series depends on the present state in first time series and previous state in second time series. Using the MCMC simulation, wind powers are generated interchangeably for both series, using matrix Px for the first site and Py for the second one. The transitions and the simulation flow are illustrated in Fig. 11. This simulation method provides significant improvement in the resulting aggregate wind power data series. In Fig. 12, it is shown that the cdf of the simulated total wind power perfectly matches the original one, with a correlation coefficient fully corresponding to the original (r ¼ 0.65). In comparison to the MCMC simulation with copula implementation, this method also obtains superior results for the pdf of 10 min variations shown in Fig. 13. Hence, unlike other simulation methods analyzed here, using the proposed MCMC approach it is possible to generate wind data almost fully corresponding to the original one, allowing satisfactory representation of short term wind power variations inside both time series, along with their aggregate power. This simulation method can generally be analogously applied to more wind farm sites, where the order of the Markov chain needs to be at least equal to number of considered sites. Transition probabilities for n matrices can be written as:
(16)
Hence, a future state in the first time series depends on the previous states in both the first and second time series.
Probability
0,7
0,6
proposed 2nd order MCMC original series
0,5
independent MCMC (ρ=0) copula generated
0,4
0,3
0,2
0,1
0,0 -10%
-5%
0%
5%
10%
Total wind power variations [% of capacity]
Fig. 11. Illustration of the proposed second-order MCMC simulation.
Fig. 13. The pdf for 10 min variations of the original and proposed 2nd order MCMC generated total wind power.
R. Goic et al. / Renewable Energy 35 (2010) 2602e2609
1st site x1-i
t
px1-ij...kl
2nd site x2-j
7. Conclusion
nth site
...
xn-k
px2-j...klm x1-l
t+1
...
x1-m
.. .
xn-o
pxn-klm...o
.. .
.. .
Fig. 14. Illustration of the proposed nth-order MCMC simulation.
P X1
tþ1
¼ x1l X1
¼ x1i ; X2
t
t
¼ x2j ; .; Xn
t
¼ xnk
¼ px1ij.kl (17) P X2
tþ1
¼ x2m X2
t
¼ x2j ; .; Xn
t
¼ xnk ; X1
tþ1
¼ x1l
¼ px2j.klm (18)
.
PfXn
tþ1
¼ xno jXn
t
¼ xnk ; X1
¼ x2m ; .g ¼ pxn
tþ1
klm.o
¼ x1m ; X2
2609
Analysis of the impact of wind power variations on power system operation requires adequate modeling of total power output from geographically dispersed wind farms. In this paper, different simulation methods are investigated on annual 10 min wind data series for two real wind sites. The obtained results show that adequate stochastic dependence modeling is required in order to obtain correct distributions of aggregate power output and total power variations. Additionally, for fitting the pdf of short-term variations, it is necessary to preserve at least short-term persistence inside generated power time series. A modified second-order Markov chain Monte Carlo simulation method is proposed in order to obtain synthetic time series of aggregate wind power which fits the original data very well with respect to the basic statistical parameters. This method can show certain limitations in considerable memory usage based on the number of considered wind sites, but does not play a major role in the preprocessing and main algorithm complexity. For future work it is planned to additionally investigate the possibility of using this method for aggregate wind farms generation modeling. Furthermore, it will be assessed whether the proposed method can also be used to model the dependency between wind power generation and load.
tþ1
(19)
The mutual dependence between wind data series at multiple sites could be simply examined using a correlation matrix containing all combinations of cross-correlations. Illustration of the simulation algorithm for n sites is shown in Fig. 14. The main limitation of the presented simulation method lies in potentially high memory space requirements. Namely, since the order of the Markov chain is equal to the number of time series, the size of the transition probability matrix exponentially grows with number of wind farms sites. Since high resolution of the states is required (as described in Subsection 4.1), the dimensions of the transition matrices can be very large, even for a relatively low number of sites. However, the complexity of the simulation approach is equal to O(m n), where m is the number of states, and n is the number of values in the data series. This complexity also refers to the preprocessing transition matrices generation step. Consequently, the number of considered sites does not play a major role in the execution time of the preprocessing and the simulation itself, although it does significantly increase the required memory usage. Furthermore, the size of the transition probability matrix, i.e. the memory usage, can significantly be reduced if short-term wind power variations are not considered the main simulation objective. It follows that the granularity of the states defined for the considered number of sites, creates a trade-off between simulation complexity and accuracy.
References [1] World Wind Energy conference and exhibition. World wind energy report 2008, www.wwindea.org/home/images/stories/worldwindenergyreport2008 _s.pdf; June 2009. [2] Papaefthymiou G, Klöckl B. MCMC for wind power simulation. IEEE Transactions on Energy Conversion 2008;23(1). [3] Ettoumi FY, Sauvageot H, Adane A-E-H. Statistical bivariate modeling of wind using first-order Markov chain and Weibull distribution. Renewable Energy 2003;28:1787e802. [4] Shamshad A, Bawadi MA, Hussin W, Majid T, Sanusi S. First and second order Markov chain models for synthetic generation of wind speed time series. Energy 2005;30:693e708. [5] Sahin AD, Sen Z. First-order Markov chain approach to wind speed modeling. Journal of Wind Engineering and Industrial Aerodynamics 2001;89:263e9. [6] Van der Hoven I. Power spectrum of horizontal wind speed in the frequency range from 0.0007 to 900 cycles per hour. Journal of Meteorology 1956;14. [7] Nørgård P, Giebel G, Holttinen H, Söder L, Petterteig A. Fluctuations and predictability of wind and hydropower, Risø-R-1443, 2004. [8] Goic R, Jakus D, Penovic I, Vasilj J, Kovacevic A, Lovric M, Klaric M. Integration of wind energy into the Croatian power system, Croatian national committee of CIGRE, Symposium on power system management, Cavtat 2008. [9] Holttinen H. Technical research Centre of Finland VTT, Impact of hourly wind Power variations on the system operation in Nordic countries, P.O. Box 1601, FIN-02044 VTT, 2005. [10] European wind energy association (EWEA). Large scale integration of wind energy in the European power supply: analysis, issues and recommendations, EWEA report, http://www.ewea.org/fileadmin/ewea_documents/documents/ publications/grid/051215_Grid_report.pdf; 2005. [11] Meyn SP, Tweedie RL. Markov chains and stochastic stability. first ed. London: Springer-Verlag; 1993. [12] Papaefthymiou G. Integration of stochastic generation in power systems, Ph.D. dissertation, Delft university of technology, ISBN: 978-90-8570-186-6; 2007. [13] Fermanian J, Radulovic D, Wegkamp M. Weak convergence of empirical copula processes. Bernoulli 2004;10(5):847e60.