- Email: [email protected]
A novel energy pattern factor method for wind speed distribution parameter estimation
Energy Conversion and Management 106 (2015) 1124–1133
Contents lists available at ScienceDirect
Energy Conversion and Management journal homepage: www.elsevier.com/locate/enconman
A novel energy pattern factor method for wind speed distribution parameter estimation Seyit Ahmet Akdag˘ ⇑, Önder Güler Istanbul Technical University, Energy Institute, Ayazaga Campus, 34469 Maslak, Istanbul, Turkey
a r t i c l e
i n f o
Article history: Received 12 June 2015 Accepted 13 October 2015 Available online 11 November 2015 Keywords: Wind energy Wind speed distribution Weibull distribution Novel energy pattern method
a b s t r a c t Power output of wind turbine depends on many factors. Among them, the most crucial one is wind speed. Since wind speed data is a significant factor for wind energy analyses, it should be modeled accurately. Weibull distribution has been used extensively to model variation of wind speed. Therefore, the most appropriate distribution parameter estimation method selection is critical in order to minimize data set modeling errors. In this context, a novel, robust, efficient and better method than standard methods to estimate Weibull parameters is presented for the first time in this paper. The accuracy of the proposed method is verified using different data sets. Also, developed method is compared with Graphic Method (GM), Maximum Likelihood Method (MLM), Alternative Maximum Likelihood Method (AMLH), Modified Maximum Likelihood Method (MMLH), Moment Method (MM), Justus Moment Method (JMM), WAsP Method (WM) and Power Density Method (PD). The results indicate that the proposed novel method is adequate to determine Weibull distribution parameters. Ó 2015 Elsevier Ltd. All rights reserved.
1. Introduction Fossil energy sources have served as the primary energy source and supplied about 86% of the world’s final energy demand [1]. According to the statistics, the current trend of consuming scarce fossil energy sources causes these energy sources deplete within the next decades [2]. Due to fluctuating fossil energy source prices, increasing consumption of limited reserves, global population rise, growth in energy demand and energy supply security problems, many countries have been forced to use renewable energy sources. Installed wind power capacity accounts for 5% of global electricity generation, equal to the 480 TW h per year by the end of 2014. Wind energy’s contribution to the global electricity supply is expected to reach 12% by 2020 and 22% by 2030 [3]. Wind energy has a technical potential to supply global electricity consumption. Wind turbine power output depends on many factors such as mean wind speed, power density, wind speed distribution, turbine hub height, turbine rated power, shape of power curve, air density, turbulence intensity, and other factors. Among them, the most crucial one is wind speed distribution. As stated by Morgan et al. [4] main uncertainty in estimation of wind turbine annual output lies in the selection of accurate distribution, since power curve of wind ⇑ Corresponding author. E-mail addresses: [email protected] (S.A. Akdag˘), [email protected] (Ö. Güler). http://dx.doi.org/10.1016/j.enconman.2015.10.042 0196-8904/Ó 2015 Elsevier Ltd. All rights reserved.
turbine is known fairly accurate. According to the International Electrotechnical Comission Standard [5] and other international recommendations, Weibull distribution is the most appropriate and widely used one to determine wind energy potential. Moreover, Weibull distribution is the default option to estimate energy output of wind turbine or wind farm for numerous wind energy softwares [6,7]. Also, Weibull distribution is recommended as a main distribution for wind analysis in wind energy textbooks. Therefore, when a new distribution is proposed to describe wind speed distribution, it is often compared with commonly used Weibull distribution. Weibull distribution is not a universal distribution to represent wind distribution for all geographical locations in the world owing to the sharp differences in climate and topography [8,9]. In recent years, entropy based distribution and mixture distributions such as normal-Weibull, Weibull–Weibull distribution have been proposed as alternatives to the Weibull distribution in literature [9–13]. The main objective of the present article is to introduce a novel, robust, efficient, practical, empirical and better method than standard methods to estimate Weibull distribution parameters. The remainder of the paper is organized as follows. In section two, literature review for parameter estimation method was carried out. In Section 3, widely used eight methods to estimate parameters of Weibull distribution were revisited. In Section 4, a novel method was introduced. In Section 5, in order to verify the suitability of this proposed novel method, the introduced method was
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1125
Nomenclature AMLH CDF GM GOF JMM MM MMLH MLM MSV NEPFM PD PDF RMSE WM k c an
Alternative Maximum Likelihood Method cumulative density function Graphic Method goodness of fit test Justus Moment Method Moment Method Modified Maximum Likelihood Method Maximum Likelihood Method Mean Sum Variation novel energy pattern factor method Power Density Method Probability density function Root mean square error WAsP Method Weibull distribution shape parameter (–) Weibull distribution scale parameter (m/s) coefficient for shape parameter
compared with eight methods through Monte Carlo Simulations. This novel method was also compared with eight methods used in terms of goodness of fit tests (GOF). Then, the influence of parameter estimation methods on the wind characteristics of the 13 selected regions in Turkey was analyzed. Section 6 concludes the paper. 2. Literature review The Weibull distribution was first introduced as Wallodi Weibull [14]. Utilization area of this distribution is vast and encompasses nearly all research areas. Thus, parameter estimation is considered as a critical topic due to the accuracy of feasibility and reliability analysis. Graphic Method (GM), Maximum Likelihood Method (MLM), Alternative Maximum Likelihood Method (AMLH), Modified Maximum Likelihood Method (MMLH), Moment Method (MM), Justus Moment Method (JMM), WAsP Method (WM) and Power Density Method (PD) are the most commonly used methods. Each of these methods has their own benefits and drawbacks. In literature, these methods are compared several times in order to investigate their efficiency, accuracy and capability in modeling data set [15–22]. It is well known that in some cases these methods can provide dramatically different results in estimation the distribution parameters. Dorvlo [16] used the chi-square, moments and regression methods to determine Weibull distribution parameters and concluded that the moment and regression methods are not as well as chi-square method. The performance of GM, JMM and MLM methods was compared by Akdag and Guler [17] considering two goodness of fit tests. It was concluded that JMM is better than other methods. PD was introduced by Akdag and Dinler [21] to estimate Weibull distribution parameters for wind energy applications. This method was compared with GM, JMM and MLM methods. Result of the study revealed that PD is adequate and able to provide high accuracy for estimation of Weibull parameters. Saleh et al. [22], discussed four methods for distribution parameter estimation, namely, MLM, MMLM, GM and PD methods. Based on the goodness of fit test results, MLM was recommended. Chang compared six methods and concluded that if wind speed distribution matches well with Weibull distribution, these methods are applicable [23]. Also, a comparison among seven methods to estimate Weibull distribution parameter was presented, as regards to their accuracy and capability to model the wind speed distributions in Brazil [24].
bn cn dn Epf Fðv Þ f ðv Þ mn
r
vm
Cð:Þ
v3
ðv m Þ3 R2 mn mdn stdn
coefficient for shape parameter coefficient for scale parameter coefficient for scale parameter energy pattern factor cumulative distribution probability density function Nth moment of Weibull distribution standard deviation mean wind speed Gamma function mean of wind speed cubes and cube of mean wind speed. coefficient of correlation Nth moment of Weibull distribution Nth moment of data sets Weibull distribution standard deviation and stdd is data set standard deviation
Table 1 summaries the goodness of fit results of the seven methods for selected regions of Brazil [24]. Methods are ranked according to three goodness of fit criteria. The rankings were done by considering maximum coefficient of determination, minimum root mean square error and minimum chi square value. GM is the best method for three cases, second best method for one case and worst for one case. PD method is the best method for two cases, second good method for three cases and third method for one case. EEM is the best method for one case, second best method for one case, third best method for two cases and worst method for one case. It was revealed that GM method is the best and followed by PD and EEM. The performance of seven parameter determination method was compared by Azad et al. [25] to find the best method in terms of six GOF tests. Results of the paper summarized in Table 2. The rankings were done by considering six GOF test results. MM method is best method for two case, second best method for one case and third method for two cases. PD is best method for two cases and second best method for two cases. It was revealed that MM method is the best and followed by MLM and PD. According to Tables 1 and 2 none of the methods ranked as the best for all cases. However, some of the methods may perform better than others. It can be concluded that there is not a single, universally accepted, best method to estimate Weibull distribution parameters. So, these literature show us that this topic is still open to exploration. In this paper, a novel method was outlined which can be efficiently and accurately used to determine Weibull parameters. Also proposed methods can be used for various distributions.
3. Methods for estimating Weibull parameter A statistical distribution to show wind characteristics provides information about the wind regime at a measurement site. This Table 1 Rank of parameter determination methods for Brazil [24].
GM PD EEM MLM MM MMLM JMM
1th
2nd
3rd
4th
5th
6th
7th
3 2 1 – – – –
1 3 1 1 – – –
– 1 2 2 1 – –
– – – 2 1 2 1
1 – – – 2 1 2
– – 1 – 1 2 2
1 – 1 1 1 1 1
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1126
Table 2 Rank of parameter determination methods Bangladesh [25].
GM PD EEM MLM MM MMLM JMM
1th
2nd
3rd
4th
5th
6th
7th
– 2 – 1 3 – –
– 1 – 3 1 1 –
– – 1 1 2 2 –
2 2 1 – – 1 –
1 1 3 – – 1 –
3 – 1 1 – 1 –
– – – – – – 6
Pn c¼
is essential both for making a decision about the economic feasibility of a wind power project and for the selection of the most suitable site for installing a wind farm through the estimation of the capacity factor. Cumulative distribution and probability density functions of Weibull distribution are given by the following Eqs. (1) and (2), respectively. v
Fðv Þ ¼ 1 eð c Þ f ðv Þ ¼
k
ð1Þ
k v k1 ðvc Þk e ; c c
v > 0; k > 0; c > 0
ð2Þ
where F(v) is cumulative Weibull distribution function, f(v) is the Weibull probability density function, k is the dimensionless shape parameter and c (in same unit of v) is the scale parameter of Weibull distribution. Abscissa of Weibull distribution is controlled with the scale parameter, in other words, scale parameter is related to the mean of data set. Shape parameter describes the flatness of Weibull distribution. If shape parameter gets larger, distribution gets narrower and its peak value gets higher. Weibull distribution can be skewed to the right side or to the left side according to the parameters. Figs. 1 and 2 show how the Weibull distribution changes for selected wide range of value of the parameters. Fig. 3 shows the power density variations with Weibull parameters. Weibull distribution is capable of approximating several distributions. For instance, Weibull distribution becomes exponential distribution when the shape parameter is equal to 1. In a special case, where k = 2 Rayleigh distribution is obtained and approximated Gaussian distribution obtained for k = 3.6. Also, When k < 1, Weibull distribution decreases monotonically and bell shaped if k > 1.
In this section, a brief description of the eight parameter estimation methods is given. 3.1.1. Maximum Likelihood Method The idea behind the MLM method is to determine the Weibull distribution parameters which maximize the likelihood function. The likelihood function for Weibull distribution is given by the following Eq. (3) n Y
f ðv i ; k; cÞ
ð3Þ
i¼1
Taking the logarithm of the likelihood function, the following Eq. (4) is obtained n X ln L ¼ nðln k k ln cÞ þ ðk 1Þ lnðv i Þ i¼1
i¼1
n 1X ðv i Þk ck i¼1
!
Maximizing Eq. (4) is achieved by computing the derivative with respect to the scale and shape parameters and setting them to zero. Distribution parameters are determined according to the equations below
1k
ð6Þ
where v i is the wind speed in time site i and n is the number of nonzero wind speed data. According to the Eq. (6), scale parameter of Weibull distribution depends on shape parameter. MLM requires significant computational effort in order to overcome this problem and so as to minimize the computational requirements of this method, AMLM and MMLM methods were proposed [26]. As stated by Kantar and Sßenog˘lu [27] If k < 1, MLM estimators do not exist because the likelihood function is not bounded, and for 1 < k < 2, MLMs exist but do not satisfy the usual regularity conditions and lead to inefficient estimators. MLM is efficient for large data size and if the data size is smaller than 250, this method is not recommended to use. According to the Chang [23], if data number becomes larger, the performance for parameter estimation methods is improved. So, this method may lose its accuracy for small data set. 3.1.2. Alternative Maximum Likelihood Method Due to the iterative characteristics of Eq. (5), MLM was simplified considering some assumptions [26]. Simple estimating procedure has been developed and called AMLH. In this method, Weibull distribution shape parameter is determined by Eq. (7)
0
10:5 nðn 1Þ P 2 A 2 n ln v ln v i i i¼1 i¼1
p k ¼ pffiffiffi @ P n 6
n
ð7Þ
Shape parameter determined with Eq. (6). 3.1.3. Modified Maximum Likelihood Method Frequency distribution format is widely used data format for statistics. So, MMLM was proposed to improve the accuracy of MLM method and to utilize frequency format. This method can be regarded as kind of weighed MLM. Shape and scale parameters are determined by following equation
Pnb
v ki ln v i f ðv i Þ Pnb k i¼1 v i f ðv i Þ
i¼1
Pnb
ln v i f ðv i Þ f ðv P 0Þ
i¼1
!1 ð8Þ
!1
Pnb c¼
i¼1
v ki f ðv i Þ
k
f ðv P 0Þ
ð9Þ
where v i is the wind data central to bin i, nb the number of bins. f ðv i Þ is the frequency for data within bin i. 3.1.4. Moment method Main idea behind the Moment method is matching data set first moment and standard deviation to the distribution moment and standard deviation. Weibull distribution nth moment is calculated by following Eq. (10).
Z mn ¼
1
0
ð4Þ
v ki
ð5Þ
n
k¼
3.1. Parameter determination methods
Lðv 1 ; v 2 . . . :v n; k; cÞ ¼
!1 Pn k Pn ln v i i¼1 v i ln v i Pn k i¼1 n i¼1 v i
k¼
v n f ðv Þdv
ð10Þ
After manipulations scale parameter of Weibull distribution can be calculated as follows
c¼
vm
C 1 þ 1k
in which CðÞ is the Gamma function.
ð11Þ
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1127
Fig. 1. Weibull distribution for shape parameter is 2 and scale parameter between 5 and 15.
Fig. 2. Weibull distribution for scale parameter is 10 and shape parameter between 1 and 15.
Fig. 3. Power density variations according to the distribution parameters.
Standard deviation shows the regularity of the distribution. Standard deviation of Weibull distribution is expressed using Eq. (12).
r¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðm2 v 2m Þ
ð12Þ
Eq. (13) shows the ratio of standard deviation to the mean value
r
vm
vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi u uC 1 þ 2 C2 1 þ 1 k k ¼t C2 1 þ 1k
3.1.5. Justus Moment Method The Justus approximation [28] is considered as the best alternative to the numerical solution of Eq. (13). According to the Justus approximation Weibull shape parameter can be determined as following equation.
k¼
r
1:086
vm
ð14Þ
Scale parameter determined with Eq. (11).
ð13Þ
Above Eq. (13) can be solved by using numerical methods to find shape parameter.
3.1.6. Graphic Method Graphic method is derived by a logarithmic function of the cumulative Weibull distribution expressed in Eq. (1). By taking the logarithm twice of Eq. (1), following equation can be obtained
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1128
Fig. 4. Epf variation for shape parameter is between 1 and 15.
Table 3 NEPFM coefficients for shape parameter. a0 a1 a2 a3 a4
0.220374 3.27527 5.78961 2.15143 0.590396
b0 b1 b2 b3 b4
1.27285 3.69115 2.60973 0.800468 0.992007
3.1.8. PD method Akdag and Dinler developed PD method considering the variation of Energy pattern factor ðEpf Þ between 1.45 and 4.4 [21]. Result of the study showed that PD method approximation works with high degree of accuracy for variation of shape parameter between 1.2 and 2.75. Main idea behind the PD Method is matching data set first and third moments to the distribution first and third moments.
k¼1þ Table 4 GOF test results for NEPFM shape parameter.
ðEpf Þ2
ð20Þ
Scale parameter determined with Eq. (11).
GOF test
Result
MEð%Þ MaxEð%Þ MinEð%Þ RMSE
0.003062 0.115717 2.07 ⁄ 106 5.62 ⁄ 109
ln½ ln½1 Fðv Þ ¼ k lnv k lnc
3:69
ð15Þ
y ¼ ax þ b line is obtained, in which a ¼ k and b ¼ klnc: Hence,
4. A novel energy pattern factor method Wind power density and its distribution are one of the most significant criteria for wind energy analysis and wind turbine design. Energy pattern factor is expressed as the ratio between mean of cubic wind speed data set and the cube of mean wind speed data set. Epf can be expressed by the following Eq. (21).
v3 ðv m Þ3
k¼a
ð16Þ
Epf ¼
b c ¼ exp a
ð17Þ
where v 3 is mean of wind speed cubes and ðv m Þ3 is cube of mean wind speed. Epf is related with the mean wind speed, power density and wind speed frequency distribution. Energy pattern factor is a sign of wind speed variability. Variation of Epf for several wind speed measurement stations and various wind regimes was reported by Earnest [29]. The value of Epf for polar regions is 6, for continental and irregular regions is 2.7, for coastal regions 1.57–1.92 and 1.22–1.36 for trade wind regions. Epf for Rayleigh distribution is calculated as 1.91. According the Weibull distribution moments Epf can be expressed with following Eq. (22),
3.1.7. WAsP Method Main idea behind the WAsP method is matching data set third moment to the distribution third moment. With the utilization of cumulative distribution function and Weibull distribution first moment following Eq. (18) is obtained.
k 1 lnð1 Fðv m ÞÞ ¼ C 1 þ k
ð18Þ
When we match the third moment of Weibull distribution to the data set third moment, following Eq. (19) obtained.
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi v 3 3 c¼ C 1 þ 3k
R1 Epf ¼ R01 0
ð19Þ
v 3 f ðv Þdv v f ðv Þdv 3
ð21Þ
ð22Þ
After some integrations and manipulations following Eq. (23) is obtained for Epf .
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1129
Fig. 5. Error (%) value of NEPFM for shape parameter.
Fig. 6. Variation of C 1 þ 1k versus k parameter between 1 and 15. Fig. 7. Errors for C 1 þ 1k versus k fitting. Table 5 GOF test results for C 1 þ 1k fitting.
Table 6 NEPFM scale parameter coefficients.
GOF test
Result
MEð%Þ MaxEð%Þ MinEð%Þ RMSE
0.03354 0.639488 5.09 ⁄ 106 3.98 ⁄ 104
k
C 1 þ 3k Epf ¼ 3 C 1 þ 1k
c0 c1
ð23Þ
Fig. 4 shows the Epf variation for shape parameter is between 1 and 15. Weibull distribution shape parameters are generally greater than 1 for wind speed distributions. Therefore, in this paper, a novel and simple approximation with higher degree of accuracy is proposed considering variation of shape parameter between 1 and 15. A formulation can be done with curve fitting technique. This variation can be modeled with Eq. (24). This method is called a novel energy pattern factor method (NEPFM). The main novelty of this novel method is its wide range of parameter estimation capability. Hence, this novel method can be used in other research areas. The expression of shape parameter according to the Epf value is given by following Eq. (24)
0.225761 0.134704
a4 E4pf þ a3 E3pf þ a2 E2pf þ a1 Epf þ a0 b4 E4pf þ b3 E3pf þ b2 E2pf þ b1 Epf þ b0
d0 d1
0.35144 0.711818
ð24Þ
where an and bn are the coefficients and they are given in Table 3. Goodness of fit test shows the suitability of proposed method. In this part of paper four GOF tests were used to show the suitability. These tests are root mean square error, absolute mean error, absolute maximum error and absolute minimum error. Also, error test show the graphical suitability. Formulation of these tests are given with following Eqs. (25–29).
MEð%Þ ¼
kn 100 X jkiestimated ki j kn i¼1 ki
MinEð%Þ ¼ 100 Min
jkiestimated ki j ki
ð25Þ
ð26Þ
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1130 Table 7 MSV (%) results for NEPFM. k (–)
c (m/s)
1.25 1.5 1.75 2 2.25 2.5 2.75 3 4 5 6 7 8 9 10
5
6
7
8
9
10
0.060945 0.079014 0.070184 0.02742 0.05076 0.010154 0.013508 0.007654 0.016538 0.010358 0.01272 0.018828 0.027866 0.018828 0.027866
0.301568 0.060089 0.031533 0.020391 0.002813 0.026924 0.005307 0.011233 0.006271 0.009126 0.012711 0.018914 0.028198 0.018914 0.028198
0.041089 0.04904 0.027042 0.055753 0.011992 0.004437 0.0095 0.017242 0.00641 0.007235 0.012266 0.019458 0.027722 0.019458 0.027722
0.055156 0.044688 0.071969 0.027907 0.002501 0.026924 0.00719 0.012865 0.005798 0.010104 0.012663 0.019448 0.028222 0.019448 0.028222
0.090931 0.048334 0.050686 0.034045 0.015123 0.009903 0.008742 0.003542 0.005302 0.008248 0.012459 0.019998 0.027294 0.019998 0.027294
0.056867 0.013762 0.032757 0.038073 0.020316 0.020316 0.006495 0.006101 0.006986 0.008039 0.0126 0.020334 0.027558 0.020334 0.027558
Table 8 Wind statistics of measurement stations. Site
Data number
Max speed
Mean
Measurement duration
Frequency (min)
Height (m)
Ankara Ayvalık Belen Biga Datça Gebze Gökçeada Karaman Kırıklareli Kocadag˘ Konya Raman Sinop
52,704 148,854 33,365 185,071 17,107 50,223 8734 8653 14,547 17,515 51,927 48,175 17,498
25.57 21 21.82 19.89 20.76 13.5 23.25 21.73 19.12 16.19 21.55 22.22 13.56
5.63777 4.33231 6.98462 4.83993 5.57689 3.28547 7.2489 4.30988 4.70042 3.66239 4.90025 4.02485 4.41037
2007–2008 2006–2011 2004–2007 2007–2009 2001–2002 2005–2006 2002 2005 2005–2007 1999–2000 2005–2006 2004–2005 1999–2000
10 60 60 60 60 10 60 60 60 60 10 10 30
30 10 30 10 10 10 10 10 10 10 10 10 10
Table 9 Shape parameters of measurement stations. Site
GM
MLM
AMLH
MMLH
JMM
NMM
PD
WM
NEPFM
Ankara Ayvalık Belen Biga Datça Gebze Gökçeada Karaman Kırıklareli Kocadag˘ Konya Raman Sinop
1.42592 1.35197 1.91868 1.74836 1.62462 1.75196 1.70603 1.42847 1.73553 1.83211 1.58155 1.55775 1.99487
1.37848 1.26979 1.93118 1.69292 1.50042 1.69665 1.69406 1.38797 1.75933 1.88474 1.51993 1.6149 1.9489
1.32184 1.36592 1.5869 1.72893 1.16936 1.68608 1.64245 1.33733 1.68612 1.8453 1.33904 1.63128 1.8056
1.37906 1.26899 1.94864 1.67816 1.59048 1.6854 1.69299 1.40732 1.76751 1.86831 1.56567 1.62471 1.95603
1.46315 1.26604 2.04752 1.71069 1.73014 1.73172 1.73867 1.44569 1.79475 1.91631 1.61311 1.6123 1.99849
1.44154 1.2504 2.025 1.68632 1.70572 1.7073 1.71424 1.42446 1.77034 1.89247 1.5893 1.58851 1.9754
1.47895 1.29289 2.12721 1.69878 1.80586 1.72618 1.74474 1.4454 1.75072 1.90296 1.60538 1.54228 2.00612
1.582130 1.190220 2.407438 1.626475 2.119664 1.685814 1.809683 1.515741 1.695878 1.968808 1.629129 1.469379 2.010226
1.479458 1.276034 2.115725 1.696711 1.799938 1.723161 1.741056 1.444781 1.746819 1.893774 1.605864 1.543467 1.994647
Table 10 Scale parameters of measurement stations. Site
GM
MLM
AMLH
MMLH
JMM
NMM
PD
WM
NEPFM
Ankara Ayvalık Belen Biga Datça Gebze Gökçeada Karaman Kırıklareli Kocadag˘ Konya Raman Sinop
5.85335 4.60534 7.07656 5.34536 5.62485 3.60385 7.96622 4.62493 5.37489 4.14747 5.31371 4.71683 4.88577
6.15167 4.67907 7.81115 5.42767 6.09425 3.67973 8.11263 4.71292 5.27042 4.12314 5.40551 4.49595 4.96151
6.07812 4.8002 7.53672 5.45656 5.76593 3.67412 8.05256 4.66213 5.21748 4.10391 5.23646 4.50839 4.88439
6.1537 4.68068 7.82591 5.41861 6.17546 3.68225 8.11683 4.73612 5.28134 4.11995 5.44904 4.50819 4.96921
6.2254 4.66476 7.88414 5.42673 6.25755 3.68667 8.13606 4.7513 5.28485 4.12836 5.46956 4.49225 4.9765
6.21269 4.65182 7.88293 5.42139 6.25183 3.68332 8.12886 4.74122 5.28112 4.12674 5.4621 4.48611 4.97535
6.23415 4.68557 7.88657 5.42418 6.27218 3.68594 8.13775 4.75117 5.27783 4.12748 5.4672 4.47268 4.97684
6.48800 4.39499 8.17000 5.30900 6.66300 3.64700 8.27900 4.89700 5.20400 4.18800 5.51000 4.34100 4.99000
6.234397 4.672786 7.886314 5.423651 6.271097 3.685478 8.136621 4.75088 5.277076 4.126778 5.46729 4.473014 4.976263
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1131
Table 11 MSV (%) results for measurement stations. Site
GM
MLM
AMLH
MMLH
JMM
NMM
PD
WM
NEPFM
Ankara Ayvalık Belen Biga Datça Gebze Gökçeada Karaman Kırıklareli Kocadag˘ Konya Raman Sinop
9.0150 9.6909 15.8735 4.9412 14.1619 5.5722 3.0255 4.4058 4.4247 2.9902 4.5251 10.9619 3.7233
4.4564 0.8605 2.6071 0.2090 9.1470 0.6678 1.1669 2.4836 0.5296 0.1287 2.6861 2.2348 0.6413
7.4641 4.0741 11.2456 0.5469 36.2439 0.8908 2.3643 4.9550 1.0182 0.6926 11.1602 2.5895 2.8381
4.4592 0.9179 2.3837 0.5166 5.8275 1.2602 1.2815 1.8235 0.4988 0.5117 1.3191 2.3152 0.6621
0.8607 0.8077 1.1874 0.5804 1.8263 0.4357 0.2717 0.3387 1.2656 0.5954 0.4725 2.2253 0.2162
1.3501 1.2977 1.4221 0.2509 2.2186 0.3760 0.6333 0.7482 0.5363 0.0215 0.4653 1.3232 0.3271
0.5228 1.5719 0.7607 0.1792 1.0070 0.2549 0.3703 0.3242 0.2702 0.2530 0.1828 0.5726 0.3913
1.4667 3.1251 1.7280 1.1917 2.6917 0.6931 0.6953 1.1898 0.5236 0.8337 0.2186 0.9254 0.1390
0.5249 0.5109 0.5339 0.1112 0.8358 0.1585 0.2550 0.2935 0.2167 0.0157 0.1945 0.5385 0.1296
MaxEð%Þ ¼ 100 Max
RMSE ¼
jkiestimated ki j ki
2 kn X jkiestimated ki j i¼1
Errorð%Þ ¼ 100
ð27Þ
ð28Þ
ki kiestimated ki ki
ð29Þ
where ME is absolute mean error, MinE absolute minimum error, MaxE is absolute maximum error, RMSE root mean square error and Error is the error. kiestimated is the estimated value of shape parameter with Eq. (24). Table 4 shows the GOF test results and Fig. 5 show the Error variation for shape parameter between 1 and 15. Proposed method is undoubtedly significant when only mean and third moment of the data set are available, which can easily be found in literature. Weibull distribution first moment can be calculated with the Eq. (10). Weibull distribution scale parameter can be expressed with Eq. (11). So, in this part C 1 þ 1k is modeled for shape param eter. Fig. 6 shows the variation of C 1 þ 1k versus k parameter between 1 and 15. The expression of C 1 þ 1k to the k is given by following Eq. (30) and Table 5 shows the GOF test results and Fig. 7 shows the errors for this fitting.
C 1þ
2 1 k þ c1 k þ c0 ; 2 k k þ d1 k þ d0 ;
for 1 6 k 6 15
v m ðk2 þ d1 k þ d0 Þ 2
k þ c1 k þ c 0
! 5 stdn stdd X 100 m m n dn þ MSVð%Þ ¼ 6 stdd mdn n¼1
ð32Þ
where mn is nt moment of Weibull distribution and mdn is nth moment of data sets. Also, stdn is Weibull distribution standard deviation and stdd is data set standard deviation. Results of the analyses are given in Table 7. The absolute mean discrepancy between the generated values and estimated wind speed values varies between 0.0025% and 0.30157%. Therefore, it can be deduced from Table 7 that proposed method fits very well.
ð30Þ 5.1. Comparison with measurement wind data
where cn and dn are the coefficients and they are given in Table 6. So, scale parameter can be determined with the following Eq. (31).
c
nificant as wind speed distribution during the wind energy potential assessment and turbine energy output analysis. Monte Carlo simulations are significant tools to generate large amount of wind speed data in order to show the accuracy and performance of the proposed method. Data sets were generated considering wind speed data. Generated data sets shape parameter varies between 1.25–10 and scale parameter varies between 5 and 10. In other words, minimum mean of data set is 4.43 m/s and maximum mean is 9.51 m/s. Power density of generated data sets varies between 67.84 W/m2 and 1826 W/ m2. Generated data sets can be downloaded from the website [31]. The accuracy and performance of the proposed method in estimating distribution parameter were evaluated on the basis of the Mean Sum Variation (MSV) test.
ð31Þ
5. Comparison of methods As stated by Celik, suitability of the distribution parameter estimation methods may depend on various factors such as number of data, shape of data distribution, format of data and selected goodness of fit criteria [30]. In this section, the discrepancies between the observed and the expected wind speed data sets moments and standard deviation estimation capability are taken as judgment criteria. Hence, moments of wind speed data sets are as sig-
In the second part of this section, parameter estimation methods reviewed in previous Sections 3 and 4 are applied to the measurement wind speed data sets to compare the performance of methods. Wind speed data used in this paper was recorded at 13 measurement stations in Turkey. Table 8 shows the wind statistics of the regions. Weibull distribution parameters of measurement sites are given in Tables 9 and 10. MSV test result for first three moment is given in Table 11. Wind speed distribution of selected regions are given in Fig. 8. According to the MSV results, proposed method gives best for all site. The discrepancy between the measured and estimated wind speed moments in terms of the MSV test varies between 0.15% and 0.83% for the proposed method. As a result, the proposed method proved to be the most accurate method to estimate Weibull distribution parameters.
1132
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
Fig. 8. Wind speed distribution of selected regions.
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
6. Conclusion A reasonably accurate knowledge of the Weibull distribution parameter at any wind energy site is critical to select optimum wind turbine and minimize energy generation cost. In this paper, a novel, robust, efficient, practical, empirical and better method, which can be used for other research areas as well as for wind energy applications, to estimate Weibull distribution parameters is presented. A comparison was made between NEPFM method and eight methods considering various data set and criteria. Comparison results indicate that the NEPFM method is slightly better fit than compared parameter estimation methods. Moreover, NEPFM has high degree of GOF test. So, proposed novel power density method to estimate Weibull distribution parameters can be applied to several research areas. Moreover, this novel method provides not only comparability in accuracy but it also offers a simpler calculation method than traditional methods. Therefore, NEPFM can be used as an alternative to estimate distribution parameters for various research areas. References [1] International Energy Agency (IEA). Key world energy statistics 2014. Paris, France: IEA; 2014. [2] International Energy Agency (IEA). World energy outlook 2014. Paris, France: IEA; 2014. [3] Global Wind Energy Council, Global Wind Energy Outlook; 2010. [4] Morgan EC, Lackner M, Vogel RM, Baise LG. Probability distributions for offshore wind speeds. Energy Convers Manage 2011;52:15–26. [5] IEC. 61400-12-1 1st ed 2005–12. Wind turbines – Part 12–1: power performance measurements of electricity producing wind turbines, International Electrotechnical Committee (IEC); 2005. [6] Wind Resource Assessment Software; 2015. Windographer, available at. [7] Mortensen NG, Landberg L, Troen I, Petersen EL. Wind Atlas analysis and application program (WAsP). Roskilde, Denmark: Riso National Laboratory; 1993. [8] Ouarda TBMJ, Charron C, Shin JY, Marpu PR, Al-Mandoos AH, Al-Tamimi MH, et al. Probability distributions of wind speed in the UAE. Energy Convers Manage 2015;93:414–34. [9] Carta JA, Ramírez P. Use of finite mixture distribution models in the analysis of wind energy in the Canarian Archipelago. Energy Convers Manage 2007;48:281–91. [10] Kantar YM, Usta I. Analysis of the upper-truncated Weibull distribution for wind speed. Energy Convers Manage 2015;96:81–8. [11] Celik AN, Makkawi A, Muneer T. Critical evaluation of wind speed frequency distribution functions. J Renew Sustain Energy 2010;2:013102.
1133
[12] Akdag SA, Bagiorgas HS, Mihalakakou G. Use of two-component Weibull mixtures in the analysis of wind speed in the Eastern Mediterranean. Appl Energy 2010;87:2566–73. [13] Ramirez P, Carta JA. The use of wind probability distributions derived from the maximum entropy principle in the analysis of wind energy. A case study. Energy Convers Manage 2005;47:2564–77. [14] Weibull W. A statistical theory of the strength of materials (No. 151). Generalstabens litografiska anstalts förlag; 1939. [15] Kantar YM, Kurban M, Hocaoglu FO. Comparison of six different parameter estimation methods in wind power applications. Sci Res Essays 2011;6:6594–604. [16] Dorvlo ASS. Estimating wind speed distribution. Energy Convers Manage 2002;43:2311–8. [17] Akdag˘ SA, Güler Ö. Calculation of wind energy potential and economic analysis by using Weibull distribution—a case study from Turkey. Part 1: determination of Weibull parameters. Energy Sources, Part B 2009;4:1–8. [18] Jowder FAL. Wind power analysis and site matching of wind turbine generators in Kingdom of Bahrain. Appl Energy 2009;86:538–45. [19] Bagiorgas HS, Giouli M, Rehman S, Al-Hadhrami LM. Weibull parameters estimation using four different methods and most energy-carrying wind speed analysis. Int J Green Energy 2011;8:529–54. [20] Maatallah T, El Alimi S, Dahmouni AW, Nasrallah SB. Wind power assessment and evaluation of electricity generation in the Gulf of Tunis, Tunisia. Sustain Cities Soc 2013;6:1–10. [21] Akdag SA, Dinler A. A new method to estimate Weibull parameters for wind energy applications. Energy Convers Manage 2009;50:1761–6. [22] Saleh H, Abou El-Azm Aly A, Abdel-Hady S. Assessment of different methods used to estimate Weibull distribution parameters for wind speed in Zafarana wind ferm, Suez Gulf. Egypt. Energy 2012;44:710–9. [23] Chang TP. Performance comparison of six numerical methods in estimating Weibull parameters for wind energy application. Appl Energy 2011;88:272–82. [24] de Andrade CF, Neto HFM, Rocha PAC, da Silva MEV. An efficiency comparison of numerical methods for determining Weibull parameters for wind energy applications: a new approach applied to the northeast region of Brazil. Energy Convers Manage 2014;86:801–8. [25] Azad AK, Rasul MG, Yusaf T. Statistical diagnosis of the best weibull methods for wind power assessment for agricultural applications. Energies 2014;7:3056–85. [26] Christofferson RD, Gillette DA. A simple estimator of the two-parameter Weibull distribution. J Clim Appl Meteorol 1987;26:323–5. [27] Kantar YM, Sßenog˘lu B. A comparative study for the location and scale parameters of the Weibull distribution with given shape parameter. Comput Geosci 2008;12:1900–9. [28] Justus CG, Hargraves WR, Mikhail A, Graber D. Methods for estimating wind speed frequency distributions. J Appl Meteorol 1978;17:350–3. [29] Earnest J. Wind power technology. Prentice-Hall of India Pvt. Ltd; 2014. [30] Celik AN. On the distributional parameters used in assessment of the suitability of wind speed probability density functions. Energy Convers Manage 2004;11–12:1735–47. [31] http://akademi.itu.edu.tr/gulerond/AkademikBilgiler/Projeler.
Contents lists available at ScienceDirect
Energy Conversion and Management journal homepage: www.elsevier.com/locate/enconman
A novel energy pattern factor method for wind speed distribution parameter estimation Seyit Ahmet Akdag˘ ⇑, Önder Güler Istanbul Technical University, Energy Institute, Ayazaga Campus, 34469 Maslak, Istanbul, Turkey
a r t i c l e
i n f o
Article history: Received 12 June 2015 Accepted 13 October 2015 Available online 11 November 2015 Keywords: Wind energy Wind speed distribution Weibull distribution Novel energy pattern method
a b s t r a c t Power output of wind turbine depends on many factors. Among them, the most crucial one is wind speed. Since wind speed data is a significant factor for wind energy analyses, it should be modeled accurately. Weibull distribution has been used extensively to model variation of wind speed. Therefore, the most appropriate distribution parameter estimation method selection is critical in order to minimize data set modeling errors. In this context, a novel, robust, efficient and better method than standard methods to estimate Weibull parameters is presented for the first time in this paper. The accuracy of the proposed method is verified using different data sets. Also, developed method is compared with Graphic Method (GM), Maximum Likelihood Method (MLM), Alternative Maximum Likelihood Method (AMLH), Modified Maximum Likelihood Method (MMLH), Moment Method (MM), Justus Moment Method (JMM), WAsP Method (WM) and Power Density Method (PD). The results indicate that the proposed novel method is adequate to determine Weibull distribution parameters. Ó 2015 Elsevier Ltd. All rights reserved.
1. Introduction Fossil energy sources have served as the primary energy source and supplied about 86% of the world’s final energy demand [1]. According to the statistics, the current trend of consuming scarce fossil energy sources causes these energy sources deplete within the next decades [2]. Due to fluctuating fossil energy source prices, increasing consumption of limited reserves, global population rise, growth in energy demand and energy supply security problems, many countries have been forced to use renewable energy sources. Installed wind power capacity accounts for 5% of global electricity generation, equal to the 480 TW h per year by the end of 2014. Wind energy’s contribution to the global electricity supply is expected to reach 12% by 2020 and 22% by 2030 [3]. Wind energy has a technical potential to supply global electricity consumption. Wind turbine power output depends on many factors such as mean wind speed, power density, wind speed distribution, turbine hub height, turbine rated power, shape of power curve, air density, turbulence intensity, and other factors. Among them, the most crucial one is wind speed distribution. As stated by Morgan et al. [4] main uncertainty in estimation of wind turbine annual output lies in the selection of accurate distribution, since power curve of wind ⇑ Corresponding author. E-mail addresses: [email protected] (S.A. Akdag˘), [email protected] (Ö. Güler). http://dx.doi.org/10.1016/j.enconman.2015.10.042 0196-8904/Ó 2015 Elsevier Ltd. All rights reserved.
turbine is known fairly accurate. According to the International Electrotechnical Comission Standard [5] and other international recommendations, Weibull distribution is the most appropriate and widely used one to determine wind energy potential. Moreover, Weibull distribution is the default option to estimate energy output of wind turbine or wind farm for numerous wind energy softwares [6,7]. Also, Weibull distribution is recommended as a main distribution for wind analysis in wind energy textbooks. Therefore, when a new distribution is proposed to describe wind speed distribution, it is often compared with commonly used Weibull distribution. Weibull distribution is not a universal distribution to represent wind distribution for all geographical locations in the world owing to the sharp differences in climate and topography [8,9]. In recent years, entropy based distribution and mixture distributions such as normal-Weibull, Weibull–Weibull distribution have been proposed as alternatives to the Weibull distribution in literature [9–13]. The main objective of the present article is to introduce a novel, robust, efficient, practical, empirical and better method than standard methods to estimate Weibull distribution parameters. The remainder of the paper is organized as follows. In section two, literature review for parameter estimation method was carried out. In Section 3, widely used eight methods to estimate parameters of Weibull distribution were revisited. In Section 4, a novel method was introduced. In Section 5, in order to verify the suitability of this proposed novel method, the introduced method was
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1125
Nomenclature AMLH CDF GM GOF JMM MM MMLH MLM MSV NEPFM PD PDF RMSE WM k c an
Alternative Maximum Likelihood Method cumulative density function Graphic Method goodness of fit test Justus Moment Method Moment Method Modified Maximum Likelihood Method Maximum Likelihood Method Mean Sum Variation novel energy pattern factor method Power Density Method Probability density function Root mean square error WAsP Method Weibull distribution shape parameter (–) Weibull distribution scale parameter (m/s) coefficient for shape parameter
compared with eight methods through Monte Carlo Simulations. This novel method was also compared with eight methods used in terms of goodness of fit tests (GOF). Then, the influence of parameter estimation methods on the wind characteristics of the 13 selected regions in Turkey was analyzed. Section 6 concludes the paper. 2. Literature review The Weibull distribution was first introduced as Wallodi Weibull [14]. Utilization area of this distribution is vast and encompasses nearly all research areas. Thus, parameter estimation is considered as a critical topic due to the accuracy of feasibility and reliability analysis. Graphic Method (GM), Maximum Likelihood Method (MLM), Alternative Maximum Likelihood Method (AMLH), Modified Maximum Likelihood Method (MMLH), Moment Method (MM), Justus Moment Method (JMM), WAsP Method (WM) and Power Density Method (PD) are the most commonly used methods. Each of these methods has their own benefits and drawbacks. In literature, these methods are compared several times in order to investigate their efficiency, accuracy and capability in modeling data set [15–22]. It is well known that in some cases these methods can provide dramatically different results in estimation the distribution parameters. Dorvlo [16] used the chi-square, moments and regression methods to determine Weibull distribution parameters and concluded that the moment and regression methods are not as well as chi-square method. The performance of GM, JMM and MLM methods was compared by Akdag and Guler [17] considering two goodness of fit tests. It was concluded that JMM is better than other methods. PD was introduced by Akdag and Dinler [21] to estimate Weibull distribution parameters for wind energy applications. This method was compared with GM, JMM and MLM methods. Result of the study revealed that PD is adequate and able to provide high accuracy for estimation of Weibull parameters. Saleh et al. [22], discussed four methods for distribution parameter estimation, namely, MLM, MMLM, GM and PD methods. Based on the goodness of fit test results, MLM was recommended. Chang compared six methods and concluded that if wind speed distribution matches well with Weibull distribution, these methods are applicable [23]. Also, a comparison among seven methods to estimate Weibull distribution parameter was presented, as regards to their accuracy and capability to model the wind speed distributions in Brazil [24].
bn cn dn Epf Fðv Þ f ðv Þ mn
r
vm
Cð:Þ
v3
ðv m Þ3 R2 mn mdn stdn
coefficient for shape parameter coefficient for scale parameter coefficient for scale parameter energy pattern factor cumulative distribution probability density function Nth moment of Weibull distribution standard deviation mean wind speed Gamma function mean of wind speed cubes and cube of mean wind speed. coefficient of correlation Nth moment of Weibull distribution Nth moment of data sets Weibull distribution standard deviation and stdd is data set standard deviation
Table 1 summaries the goodness of fit results of the seven methods for selected regions of Brazil [24]. Methods are ranked according to three goodness of fit criteria. The rankings were done by considering maximum coefficient of determination, minimum root mean square error and minimum chi square value. GM is the best method for three cases, second best method for one case and worst for one case. PD method is the best method for two cases, second good method for three cases and third method for one case. EEM is the best method for one case, second best method for one case, third best method for two cases and worst method for one case. It was revealed that GM method is the best and followed by PD and EEM. The performance of seven parameter determination method was compared by Azad et al. [25] to find the best method in terms of six GOF tests. Results of the paper summarized in Table 2. The rankings were done by considering six GOF test results. MM method is best method for two case, second best method for one case and third method for two cases. PD is best method for two cases and second best method for two cases. It was revealed that MM method is the best and followed by MLM and PD. According to Tables 1 and 2 none of the methods ranked as the best for all cases. However, some of the methods may perform better than others. It can be concluded that there is not a single, universally accepted, best method to estimate Weibull distribution parameters. So, these literature show us that this topic is still open to exploration. In this paper, a novel method was outlined which can be efficiently and accurately used to determine Weibull parameters. Also proposed methods can be used for various distributions.
3. Methods for estimating Weibull parameter A statistical distribution to show wind characteristics provides information about the wind regime at a measurement site. This Table 1 Rank of parameter determination methods for Brazil [24].
GM PD EEM MLM MM MMLM JMM
1th
2nd
3rd
4th
5th
6th
7th
3 2 1 – – – –
1 3 1 1 – – –
– 1 2 2 1 – –
– – – 2 1 2 1
1 – – – 2 1 2
– – 1 – 1 2 2
1 – 1 1 1 1 1
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1126
Table 2 Rank of parameter determination methods Bangladesh [25].
GM PD EEM MLM MM MMLM JMM
1th
2nd
3rd
4th
5th
6th
7th
– 2 – 1 3 – –
– 1 – 3 1 1 –
– – 1 1 2 2 –
2 2 1 – – 1 –
1 1 3 – – 1 –
3 – 1 1 – 1 –
– – – – – – 6
Pn c¼
is essential both for making a decision about the economic feasibility of a wind power project and for the selection of the most suitable site for installing a wind farm through the estimation of the capacity factor. Cumulative distribution and probability density functions of Weibull distribution are given by the following Eqs. (1) and (2), respectively. v
Fðv Þ ¼ 1 eð c Þ f ðv Þ ¼
k
ð1Þ
k v k1 ðvc Þk e ; c c
v > 0; k > 0; c > 0
ð2Þ
where F(v) is cumulative Weibull distribution function, f(v) is the Weibull probability density function, k is the dimensionless shape parameter and c (in same unit of v) is the scale parameter of Weibull distribution. Abscissa of Weibull distribution is controlled with the scale parameter, in other words, scale parameter is related to the mean of data set. Shape parameter describes the flatness of Weibull distribution. If shape parameter gets larger, distribution gets narrower and its peak value gets higher. Weibull distribution can be skewed to the right side or to the left side according to the parameters. Figs. 1 and 2 show how the Weibull distribution changes for selected wide range of value of the parameters. Fig. 3 shows the power density variations with Weibull parameters. Weibull distribution is capable of approximating several distributions. For instance, Weibull distribution becomes exponential distribution when the shape parameter is equal to 1. In a special case, where k = 2 Rayleigh distribution is obtained and approximated Gaussian distribution obtained for k = 3.6. Also, When k < 1, Weibull distribution decreases monotonically and bell shaped if k > 1.
In this section, a brief description of the eight parameter estimation methods is given. 3.1.1. Maximum Likelihood Method The idea behind the MLM method is to determine the Weibull distribution parameters which maximize the likelihood function. The likelihood function for Weibull distribution is given by the following Eq. (3) n Y
f ðv i ; k; cÞ
ð3Þ
i¼1
Taking the logarithm of the likelihood function, the following Eq. (4) is obtained n X ln L ¼ nðln k k ln cÞ þ ðk 1Þ lnðv i Þ i¼1
i¼1
n 1X ðv i Þk ck i¼1
!
Maximizing Eq. (4) is achieved by computing the derivative with respect to the scale and shape parameters and setting them to zero. Distribution parameters are determined according to the equations below
1k
ð6Þ
where v i is the wind speed in time site i and n is the number of nonzero wind speed data. According to the Eq. (6), scale parameter of Weibull distribution depends on shape parameter. MLM requires significant computational effort in order to overcome this problem and so as to minimize the computational requirements of this method, AMLM and MMLM methods were proposed [26]. As stated by Kantar and Sßenog˘lu [27] If k < 1, MLM estimators do not exist because the likelihood function is not bounded, and for 1 < k < 2, MLMs exist but do not satisfy the usual regularity conditions and lead to inefficient estimators. MLM is efficient for large data size and if the data size is smaller than 250, this method is not recommended to use. According to the Chang [23], if data number becomes larger, the performance for parameter estimation methods is improved. So, this method may lose its accuracy for small data set. 3.1.2. Alternative Maximum Likelihood Method Due to the iterative characteristics of Eq. (5), MLM was simplified considering some assumptions [26]. Simple estimating procedure has been developed and called AMLH. In this method, Weibull distribution shape parameter is determined by Eq. (7)
0
10:5 nðn 1Þ P 2 A 2 n ln v ln v i i i¼1 i¼1
p k ¼ pffiffiffi @ P n 6
n
ð7Þ
Shape parameter determined with Eq. (6). 3.1.3. Modified Maximum Likelihood Method Frequency distribution format is widely used data format for statistics. So, MMLM was proposed to improve the accuracy of MLM method and to utilize frequency format. This method can be regarded as kind of weighed MLM. Shape and scale parameters are determined by following equation
Pnb
v ki ln v i f ðv i Þ Pnb k i¼1 v i f ðv i Þ
i¼1
Pnb
ln v i f ðv i Þ f ðv P 0Þ
i¼1
!1 ð8Þ
!1
Pnb c¼
i¼1
v ki f ðv i Þ
k
f ðv P 0Þ
ð9Þ
where v i is the wind data central to bin i, nb the number of bins. f ðv i Þ is the frequency for data within bin i. 3.1.4. Moment method Main idea behind the Moment method is matching data set first moment and standard deviation to the distribution moment and standard deviation. Weibull distribution nth moment is calculated by following Eq. (10).
Z mn ¼
1
0
ð4Þ
v ki
ð5Þ
n
k¼
3.1. Parameter determination methods
Lðv 1 ; v 2 . . . :v n; k; cÞ ¼
!1 Pn k Pn ln v i i¼1 v i ln v i Pn k i¼1 n i¼1 v i
k¼
v n f ðv Þdv
ð10Þ
After manipulations scale parameter of Weibull distribution can be calculated as follows
c¼
vm
C 1 þ 1k
in which CðÞ is the Gamma function.
ð11Þ
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1127
Fig. 1. Weibull distribution for shape parameter is 2 and scale parameter between 5 and 15.
Fig. 2. Weibull distribution for scale parameter is 10 and shape parameter between 1 and 15.
Fig. 3. Power density variations according to the distribution parameters.
Standard deviation shows the regularity of the distribution. Standard deviation of Weibull distribution is expressed using Eq. (12).
r¼
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðm2 v 2m Þ
ð12Þ
Eq. (13) shows the ratio of standard deviation to the mean value
r
vm
vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi u uC 1 þ 2 C2 1 þ 1 k k ¼t C2 1 þ 1k
3.1.5. Justus Moment Method The Justus approximation [28] is considered as the best alternative to the numerical solution of Eq. (13). According to the Justus approximation Weibull shape parameter can be determined as following equation.
k¼
r
1:086
vm
ð14Þ
Scale parameter determined with Eq. (11).
ð13Þ
Above Eq. (13) can be solved by using numerical methods to find shape parameter.
3.1.6. Graphic Method Graphic method is derived by a logarithmic function of the cumulative Weibull distribution expressed in Eq. (1). By taking the logarithm twice of Eq. (1), following equation can be obtained
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1128
Fig. 4. Epf variation for shape parameter is between 1 and 15.
Table 3 NEPFM coefficients for shape parameter. a0 a1 a2 a3 a4
0.220374 3.27527 5.78961 2.15143 0.590396
b0 b1 b2 b3 b4
1.27285 3.69115 2.60973 0.800468 0.992007
3.1.8. PD method Akdag and Dinler developed PD method considering the variation of Energy pattern factor ðEpf Þ between 1.45 and 4.4 [21]. Result of the study showed that PD method approximation works with high degree of accuracy for variation of shape parameter between 1.2 and 2.75. Main idea behind the PD Method is matching data set first and third moments to the distribution first and third moments.
k¼1þ Table 4 GOF test results for NEPFM shape parameter.
ðEpf Þ2
ð20Þ
Scale parameter determined with Eq. (11).
GOF test
Result
MEð%Þ MaxEð%Þ MinEð%Þ RMSE
0.003062 0.115717 2.07 ⁄ 106 5.62 ⁄ 109
ln½ ln½1 Fðv Þ ¼ k lnv k lnc
3:69
ð15Þ
y ¼ ax þ b line is obtained, in which a ¼ k and b ¼ klnc: Hence,
4. A novel energy pattern factor method Wind power density and its distribution are one of the most significant criteria for wind energy analysis and wind turbine design. Energy pattern factor is expressed as the ratio between mean of cubic wind speed data set and the cube of mean wind speed data set. Epf can be expressed by the following Eq. (21).
v3 ðv m Þ3
k¼a
ð16Þ
Epf ¼
b c ¼ exp a
ð17Þ
where v 3 is mean of wind speed cubes and ðv m Þ3 is cube of mean wind speed. Epf is related with the mean wind speed, power density and wind speed frequency distribution. Energy pattern factor is a sign of wind speed variability. Variation of Epf for several wind speed measurement stations and various wind regimes was reported by Earnest [29]. The value of Epf for polar regions is 6, for continental and irregular regions is 2.7, for coastal regions 1.57–1.92 and 1.22–1.36 for trade wind regions. Epf for Rayleigh distribution is calculated as 1.91. According the Weibull distribution moments Epf can be expressed with following Eq. (22),
3.1.7. WAsP Method Main idea behind the WAsP method is matching data set third moment to the distribution third moment. With the utilization of cumulative distribution function and Weibull distribution first moment following Eq. (18) is obtained.
k 1 lnð1 Fðv m ÞÞ ¼ C 1 þ k
ð18Þ
When we match the third moment of Weibull distribution to the data set third moment, following Eq. (19) obtained.
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi v 3 3 c¼ C 1 þ 3k
R1 Epf ¼ R01 0
ð19Þ
v 3 f ðv Þdv v f ðv Þdv 3
ð21Þ
ð22Þ
After some integrations and manipulations following Eq. (23) is obtained for Epf .
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1129
Fig. 5. Error (%) value of NEPFM for shape parameter.
Fig. 6. Variation of C 1 þ 1k versus k parameter between 1 and 15. Fig. 7. Errors for C 1 þ 1k versus k fitting. Table 5 GOF test results for C 1 þ 1k fitting.
Table 6 NEPFM scale parameter coefficients.
GOF test
Result
MEð%Þ MaxEð%Þ MinEð%Þ RMSE
0.03354 0.639488 5.09 ⁄ 106 3.98 ⁄ 104
k
C 1 þ 3k Epf ¼ 3 C 1 þ 1k
c0 c1
ð23Þ
Fig. 4 shows the Epf variation for shape parameter is between 1 and 15. Weibull distribution shape parameters are generally greater than 1 for wind speed distributions. Therefore, in this paper, a novel and simple approximation with higher degree of accuracy is proposed considering variation of shape parameter between 1 and 15. A formulation can be done with curve fitting technique. This variation can be modeled with Eq. (24). This method is called a novel energy pattern factor method (NEPFM). The main novelty of this novel method is its wide range of parameter estimation capability. Hence, this novel method can be used in other research areas. The expression of shape parameter according to the Epf value is given by following Eq. (24)
0.225761 0.134704
a4 E4pf þ a3 E3pf þ a2 E2pf þ a1 Epf þ a0 b4 E4pf þ b3 E3pf þ b2 E2pf þ b1 Epf þ b0
d0 d1
0.35144 0.711818
ð24Þ
where an and bn are the coefficients and they are given in Table 3. Goodness of fit test shows the suitability of proposed method. In this part of paper four GOF tests were used to show the suitability. These tests are root mean square error, absolute mean error, absolute maximum error and absolute minimum error. Also, error test show the graphical suitability. Formulation of these tests are given with following Eqs. (25–29).
MEð%Þ ¼
kn 100 X jkiestimated ki j kn i¼1 ki
MinEð%Þ ¼ 100 Min
jkiestimated ki j ki
ð25Þ
ð26Þ
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1130 Table 7 MSV (%) results for NEPFM. k (–)
c (m/s)
1.25 1.5 1.75 2 2.25 2.5 2.75 3 4 5 6 7 8 9 10
5
6
7
8
9
10
0.060945 0.079014 0.070184 0.02742 0.05076 0.010154 0.013508 0.007654 0.016538 0.010358 0.01272 0.018828 0.027866 0.018828 0.027866
0.301568 0.060089 0.031533 0.020391 0.002813 0.026924 0.005307 0.011233 0.006271 0.009126 0.012711 0.018914 0.028198 0.018914 0.028198
0.041089 0.04904 0.027042 0.055753 0.011992 0.004437 0.0095 0.017242 0.00641 0.007235 0.012266 0.019458 0.027722 0.019458 0.027722
0.055156 0.044688 0.071969 0.027907 0.002501 0.026924 0.00719 0.012865 0.005798 0.010104 0.012663 0.019448 0.028222 0.019448 0.028222
0.090931 0.048334 0.050686 0.034045 0.015123 0.009903 0.008742 0.003542 0.005302 0.008248 0.012459 0.019998 0.027294 0.019998 0.027294
0.056867 0.013762 0.032757 0.038073 0.020316 0.020316 0.006495 0.006101 0.006986 0.008039 0.0126 0.020334 0.027558 0.020334 0.027558
Table 8 Wind statistics of measurement stations. Site
Data number
Max speed
Mean
Measurement duration
Frequency (min)
Height (m)
Ankara Ayvalık Belen Biga Datça Gebze Gökçeada Karaman Kırıklareli Kocadag˘ Konya Raman Sinop
52,704 148,854 33,365 185,071 17,107 50,223 8734 8653 14,547 17,515 51,927 48,175 17,498
25.57 21 21.82 19.89 20.76 13.5 23.25 21.73 19.12 16.19 21.55 22.22 13.56
5.63777 4.33231 6.98462 4.83993 5.57689 3.28547 7.2489 4.30988 4.70042 3.66239 4.90025 4.02485 4.41037
2007–2008 2006–2011 2004–2007 2007–2009 2001–2002 2005–2006 2002 2005 2005–2007 1999–2000 2005–2006 2004–2005 1999–2000
10 60 60 60 60 10 60 60 60 60 10 10 30
30 10 30 10 10 10 10 10 10 10 10 10 10
Table 9 Shape parameters of measurement stations. Site
GM
MLM
AMLH
MMLH
JMM
NMM
PD
WM
NEPFM
Ankara Ayvalık Belen Biga Datça Gebze Gökçeada Karaman Kırıklareli Kocadag˘ Konya Raman Sinop
1.42592 1.35197 1.91868 1.74836 1.62462 1.75196 1.70603 1.42847 1.73553 1.83211 1.58155 1.55775 1.99487
1.37848 1.26979 1.93118 1.69292 1.50042 1.69665 1.69406 1.38797 1.75933 1.88474 1.51993 1.6149 1.9489
1.32184 1.36592 1.5869 1.72893 1.16936 1.68608 1.64245 1.33733 1.68612 1.8453 1.33904 1.63128 1.8056
1.37906 1.26899 1.94864 1.67816 1.59048 1.6854 1.69299 1.40732 1.76751 1.86831 1.56567 1.62471 1.95603
1.46315 1.26604 2.04752 1.71069 1.73014 1.73172 1.73867 1.44569 1.79475 1.91631 1.61311 1.6123 1.99849
1.44154 1.2504 2.025 1.68632 1.70572 1.7073 1.71424 1.42446 1.77034 1.89247 1.5893 1.58851 1.9754
1.47895 1.29289 2.12721 1.69878 1.80586 1.72618 1.74474 1.4454 1.75072 1.90296 1.60538 1.54228 2.00612
1.582130 1.190220 2.407438 1.626475 2.119664 1.685814 1.809683 1.515741 1.695878 1.968808 1.629129 1.469379 2.010226
1.479458 1.276034 2.115725 1.696711 1.799938 1.723161 1.741056 1.444781 1.746819 1.893774 1.605864 1.543467 1.994647
Table 10 Scale parameters of measurement stations. Site
GM
MLM
AMLH
MMLH
JMM
NMM
PD
WM
NEPFM
Ankara Ayvalık Belen Biga Datça Gebze Gökçeada Karaman Kırıklareli Kocadag˘ Konya Raman Sinop
5.85335 4.60534 7.07656 5.34536 5.62485 3.60385 7.96622 4.62493 5.37489 4.14747 5.31371 4.71683 4.88577
6.15167 4.67907 7.81115 5.42767 6.09425 3.67973 8.11263 4.71292 5.27042 4.12314 5.40551 4.49595 4.96151
6.07812 4.8002 7.53672 5.45656 5.76593 3.67412 8.05256 4.66213 5.21748 4.10391 5.23646 4.50839 4.88439
6.1537 4.68068 7.82591 5.41861 6.17546 3.68225 8.11683 4.73612 5.28134 4.11995 5.44904 4.50819 4.96921
6.2254 4.66476 7.88414 5.42673 6.25755 3.68667 8.13606 4.7513 5.28485 4.12836 5.46956 4.49225 4.9765
6.21269 4.65182 7.88293 5.42139 6.25183 3.68332 8.12886 4.74122 5.28112 4.12674 5.4621 4.48611 4.97535
6.23415 4.68557 7.88657 5.42418 6.27218 3.68594 8.13775 4.75117 5.27783 4.12748 5.4672 4.47268 4.97684
6.48800 4.39499 8.17000 5.30900 6.66300 3.64700 8.27900 4.89700 5.20400 4.18800 5.51000 4.34100 4.99000
6.234397 4.672786 7.886314 5.423651 6.271097 3.685478 8.136621 4.75088 5.277076 4.126778 5.46729 4.473014 4.976263
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
1131
Table 11 MSV (%) results for measurement stations. Site
GM
MLM
AMLH
MMLH
JMM
NMM
PD
WM
NEPFM
Ankara Ayvalık Belen Biga Datça Gebze Gökçeada Karaman Kırıklareli Kocadag˘ Konya Raman Sinop
9.0150 9.6909 15.8735 4.9412 14.1619 5.5722 3.0255 4.4058 4.4247 2.9902 4.5251 10.9619 3.7233
4.4564 0.8605 2.6071 0.2090 9.1470 0.6678 1.1669 2.4836 0.5296 0.1287 2.6861 2.2348 0.6413
7.4641 4.0741 11.2456 0.5469 36.2439 0.8908 2.3643 4.9550 1.0182 0.6926 11.1602 2.5895 2.8381
4.4592 0.9179 2.3837 0.5166 5.8275 1.2602 1.2815 1.8235 0.4988 0.5117 1.3191 2.3152 0.6621
0.8607 0.8077 1.1874 0.5804 1.8263 0.4357 0.2717 0.3387 1.2656 0.5954 0.4725 2.2253 0.2162
1.3501 1.2977 1.4221 0.2509 2.2186 0.3760 0.6333 0.7482 0.5363 0.0215 0.4653 1.3232 0.3271
0.5228 1.5719 0.7607 0.1792 1.0070 0.2549 0.3703 0.3242 0.2702 0.2530 0.1828 0.5726 0.3913
1.4667 3.1251 1.7280 1.1917 2.6917 0.6931 0.6953 1.1898 0.5236 0.8337 0.2186 0.9254 0.1390
0.5249 0.5109 0.5339 0.1112 0.8358 0.1585 0.2550 0.2935 0.2167 0.0157 0.1945 0.5385 0.1296
MaxEð%Þ ¼ 100 Max
RMSE ¼
jkiestimated ki j ki
2 kn X jkiestimated ki j i¼1
Errorð%Þ ¼ 100
ð27Þ
ð28Þ
ki kiestimated ki ki
ð29Þ
where ME is absolute mean error, MinE absolute minimum error, MaxE is absolute maximum error, RMSE root mean square error and Error is the error. kiestimated is the estimated value of shape parameter with Eq. (24). Table 4 shows the GOF test results and Fig. 5 show the Error variation for shape parameter between 1 and 15. Proposed method is undoubtedly significant when only mean and third moment of the data set are available, which can easily be found in literature. Weibull distribution first moment can be calculated with the Eq. (10). Weibull distribution scale parameter can be expressed with Eq. (11). So, in this part C 1 þ 1k is modeled for shape param eter. Fig. 6 shows the variation of C 1 þ 1k versus k parameter between 1 and 15. The expression of C 1 þ 1k to the k is given by following Eq. (30) and Table 5 shows the GOF test results and Fig. 7 shows the errors for this fitting.
C 1þ
2 1 k þ c1 k þ c0 ; 2 k k þ d1 k þ d0 ;
for 1 6 k 6 15
v m ðk2 þ d1 k þ d0 Þ 2
k þ c1 k þ c 0
! 5 stdn stdd X 100 m m n dn þ MSVð%Þ ¼ 6 stdd mdn n¼1
ð32Þ
where mn is nt moment of Weibull distribution and mdn is nth moment of data sets. Also, stdn is Weibull distribution standard deviation and stdd is data set standard deviation. Results of the analyses are given in Table 7. The absolute mean discrepancy between the generated values and estimated wind speed values varies between 0.0025% and 0.30157%. Therefore, it can be deduced from Table 7 that proposed method fits very well.
ð30Þ 5.1. Comparison with measurement wind data
where cn and dn are the coefficients and they are given in Table 6. So, scale parameter can be determined with the following Eq. (31).
c
nificant as wind speed distribution during the wind energy potential assessment and turbine energy output analysis. Monte Carlo simulations are significant tools to generate large amount of wind speed data in order to show the accuracy and performance of the proposed method. Data sets were generated considering wind speed data. Generated data sets shape parameter varies between 1.25–10 and scale parameter varies between 5 and 10. In other words, minimum mean of data set is 4.43 m/s and maximum mean is 9.51 m/s. Power density of generated data sets varies between 67.84 W/m2 and 1826 W/ m2. Generated data sets can be downloaded from the website [31]. The accuracy and performance of the proposed method in estimating distribution parameter were evaluated on the basis of the Mean Sum Variation (MSV) test.
ð31Þ
5. Comparison of methods As stated by Celik, suitability of the distribution parameter estimation methods may depend on various factors such as number of data, shape of data distribution, format of data and selected goodness of fit criteria [30]. In this section, the discrepancies between the observed and the expected wind speed data sets moments and standard deviation estimation capability are taken as judgment criteria. Hence, moments of wind speed data sets are as sig-
In the second part of this section, parameter estimation methods reviewed in previous Sections 3 and 4 are applied to the measurement wind speed data sets to compare the performance of methods. Wind speed data used in this paper was recorded at 13 measurement stations in Turkey. Table 8 shows the wind statistics of the regions. Weibull distribution parameters of measurement sites are given in Tables 9 and 10. MSV test result for first three moment is given in Table 11. Wind speed distribution of selected regions are given in Fig. 8. According to the MSV results, proposed method gives best for all site. The discrepancy between the measured and estimated wind speed moments in terms of the MSV test varies between 0.15% and 0.83% for the proposed method. As a result, the proposed method proved to be the most accurate method to estimate Weibull distribution parameters.
1132
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
Fig. 8. Wind speed distribution of selected regions.
S.A. Akdag˘, Ö. Güler / Energy Conversion and Management 106 (2015) 1124–1133
6. Conclusion A reasonably accurate knowledge of the Weibull distribution parameter at any wind energy site is critical to select optimum wind turbine and minimize energy generation cost. In this paper, a novel, robust, efficient, practical, empirical and better method, which can be used for other research areas as well as for wind energy applications, to estimate Weibull distribution parameters is presented. A comparison was made between NEPFM method and eight methods considering various data set and criteria. Comparison results indicate that the NEPFM method is slightly better fit than compared parameter estimation methods. Moreover, NEPFM has high degree of GOF test. So, proposed novel power density method to estimate Weibull distribution parameters can be applied to several research areas. Moreover, this novel method provides not only comparability in accuracy but it also offers a simpler calculation method than traditional methods. Therefore, NEPFM can be used as an alternative to estimate distribution parameters for various research areas. References [1] International Energy Agency (IEA). Key world energy statistics 2014. Paris, France: IEA; 2014. [2] International Energy Agency (IEA). World energy outlook 2014. Paris, France: IEA; 2014. [3] Global Wind Energy Council, Global Wind Energy Outlook; 2010. [4] Morgan EC, Lackner M, Vogel RM, Baise LG. Probability distributions for offshore wind speeds. Energy Convers Manage 2011;52:15–26. [5] IEC. 61400-12-1 1st ed 2005–12. Wind turbines – Part 12–1: power performance measurements of electricity producing wind turbines, International Electrotechnical Committee (IEC); 2005. [6] Wind Resource Assessment Software; 2015. Windographer, available at
1133
[12] Akdag SA, Bagiorgas HS, Mihalakakou G. Use of two-component Weibull mixtures in the analysis of wind speed in the Eastern Mediterranean. Appl Energy 2010;87:2566–73. [13] Ramirez P, Carta JA. The use of wind probability distributions derived from the maximum entropy principle in the analysis of wind energy. A case study. Energy Convers Manage 2005;47:2564–77. [14] Weibull W. A statistical theory of the strength of materials (No. 151). Generalstabens litografiska anstalts förlag; 1939. [15] Kantar YM, Kurban M, Hocaoglu FO. Comparison of six different parameter estimation methods in wind power applications. Sci Res Essays 2011;6:6594–604. [16] Dorvlo ASS. Estimating wind speed distribution. Energy Convers Manage 2002;43:2311–8. [17] Akdag˘ SA, Güler Ö. Calculation of wind energy potential and economic analysis by using Weibull distribution—a case study from Turkey. Part 1: determination of Weibull parameters. Energy Sources, Part B 2009;4:1–8. [18] Jowder FAL. Wind power analysis and site matching of wind turbine generators in Kingdom of Bahrain. Appl Energy 2009;86:538–45. [19] Bagiorgas HS, Giouli M, Rehman S, Al-Hadhrami LM. Weibull parameters estimation using four different methods and most energy-carrying wind speed analysis. Int J Green Energy 2011;8:529–54. [20] Maatallah T, El Alimi S, Dahmouni AW, Nasrallah SB. Wind power assessment and evaluation of electricity generation in the Gulf of Tunis, Tunisia. Sustain Cities Soc 2013;6:1–10. [21] Akdag SA, Dinler A. A new method to estimate Weibull parameters for wind energy applications. Energy Convers Manage 2009;50:1761–6. [22] Saleh H, Abou El-Azm Aly A, Abdel-Hady S. Assessment of different methods used to estimate Weibull distribution parameters for wind speed in Zafarana wind ferm, Suez Gulf. Egypt. Energy 2012;44:710–9. [23] Chang TP. Performance comparison of six numerical methods in estimating Weibull parameters for wind energy application. Appl Energy 2011;88:272–82. [24] de Andrade CF, Neto HFM, Rocha PAC, da Silva MEV. An efficiency comparison of numerical methods for determining Weibull parameters for wind energy applications: a new approach applied to the northeast region of Brazil. Energy Convers Manage 2014;86:801–8. [25] Azad AK, Rasul MG, Yusaf T. Statistical diagnosis of the best weibull methods for wind power assessment for agricultural applications. Energies 2014;7:3056–85. [26] Christofferson RD, Gillette DA. A simple estimator of the two-parameter Weibull distribution. J Clim Appl Meteorol 1987;26:323–5. [27] Kantar YM, Sßenog˘lu B. A comparative study for the location and scale parameters of the Weibull distribution with given shape parameter. Comput Geosci 2008;12:1900–9. [28] Justus CG, Hargraves WR, Mikhail A, Graber D. Methods for estimating wind speed frequency distributions. J Appl Meteorol 1978;17:350–3. [29] Earnest J. Wind power technology. Prentice-Hall of India Pvt. Ltd; 2014. [30] Celik AN. On the distributional parameters used in assessment of the suitability of wind speed probability density functions. Energy Convers Manage 2004;11–12:1735–47. [31] http://akademi.itu.edu.tr/gulerond/AkademikBilgiler/Projeler.