Extreme wind speeds from long-term synthetic records

J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

Contents lists available at SciVerse ScienceDirect

Journal of Wind Engineering and Industrial Aerodynamics journal homepage: www.elsevier.com/locate/jweia

Extreme wind speeds from long-term synthetic records Alessio Torrielli, Maria Pia Repetto n, Giovanni Solari Department of Civil, Environmental and Architectural Engineering, University of Genoa, Via Montallegro, 1, 16145 Genoa, Italy

a r t i c l e i n f o

a b s t r a c t

Article history: Received 15 June 2011 Received in revised form 11 December 2012 Accepted 22 December 2012 Available online 26 February 2013

This study proposes a novel methodology to create a large sized synthetic dataset of wind velocities and adopts this to discuss the probability distributions commonly used for extreme winds. A large number of long-term time series of mean wind speed are generated by a numerical procedure that faithfully reproduces the macro-meteorological component of wind velocity, while guaranteeing sample functions with random extremes. Through application of this technique, a large sized dataset of synthetic extreme wind observations has been extracted, of a size unprecedented in literature. Commonly applied extreme value (EV) methods are then used to process the dataset produced. In the first instance, the effectiveness of these models is tested to exclude any false effects due to the limited period covered by current wind measurements. Following this, interval estimations of design wind speeds are derived by analyzing EVs from records of different lengths in order to explore the applicability of EV distributions to real situations. The comparison between analytical and numerical results provides many interesting and intriguing points of discussion, and opens the way to new research horizons in EV analysis. & 2013 Elsevier Ltd. All rights reserved.

Keywords: Extreme value distribution Extreme wind speeds Long-term time series Monte Carlo simulation Parent distribution Wind velocity power spectrum

1. Introduction Modern wind engineering design requires the probabilistic characterization of extreme values (EV), usually, the annual maximum wind speed distribution. EV analysis is still a controversial matter, where the difference in opinions obstructs agreement on clear guidelines for practical applications. This subject became a discussion topic with the publication of the pioneering paper by Fisher and Tippett (1928) on asymptotic analysis, followed by contributions from Von Mises (1936) and Cramer (1946). Gumbel (1958) defined a clear formulation for asymptotic analysis, introducing Type I, II and III asymptotes, that encourages widely its applications. Unfortunately, no guidelines were provided on which asymptote to select in practical applications. The lack of such guidelines has led to a heated ongoing debate, which has nevertheless produced the development of alternative techniques addressing the same issue of extreme values. One of the main weaknesses of classical asymptotic analysis is the limited number of extremes involved in the analysis considering only the annual maximum values. As a result, based on the wind records currently available, the number of data items analyzed does not exceed 40–60. Poor datasets are strongly impacted by sampling errors, leaving space for serious doubts

n

Corresponding author. E-mail addresses: [email protected] (A. Torrielli), [email protected] (M.P. Repetto), [email protected] (G. Solari). 0167-6105/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.jweia.2012.12.008

on the description of the annual maximum speed distribution and resulting predictions for design values. In order to improve the accuracy of knowledge about this distribution, alternative techniques have been developed, such as process or level-crossing analysis (Gomes and Vickery, 1977), which describes the distribution of extreme values by analyzing the population of data. In parallel, other approaches were developed, which considered more than one single value per year. The peak over threshold (POT) approach (Weiss, 1971) is probably the oldest, although it has come to be applied widely only over the last two decades (Simiu and Heckert, 1996; Holmes and Moriarty, 1999; Holmes 2003). This method analyzes all values exceeding a predefined threshold. The r-LOS method (Weissman, 1978) selects the r largest observations per epoch. The method of independent storms (MIS), originally developed by Cook (1982), studies wind storm maxima, which are more frequent in number than the annual maxima. However, as observed by Palutikof et al. (1999), techniques such as POT, r-LOS and MIS require subjective decisions to be taken during the calculation; in particular, the minimum separation time between extremes, the size of r and the left-censor speed value all strongly influence the estimate of the parameters of the distribution. More recently, Cook and Harris (2004, 2008) proposed the use of FT1 penultimate distribution. This method derives from classical asymptotic analysis, taking into account finite values for the annual rate of independent events. The distribution assumes an elegant shape when the extremes are drawn from a tail equivalent Weibull parent.

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

This paper starts with an overview of the most commonly applied EV techniques (Section 2), to point out their theoretical bases, to highlight their shortcomings and to provide some guidelines on their applications. It is apparent that, notwithstanding the hot debate on the most efficient EV distribution of annual maximum wind speeds, the scientific community agrees that real observations are not of sufficient length to guarantee a satisfactory probabilistic characterization of extreme winds. Based upon these considerations, a novel approach is introduced with the aim of overcoming this problem through the generation of a great number of long-term synthetic wind speed records, simulated by means of a Monte Carlo-based procedure (Section 3). In order to apply the Monte Carlo simulation procedure (Section 4), two different blocks, each comprising three anemometer measurements, recorded in the central area of Italy, are examined. After being rendered homogeneous, these wind measurements are analyzed in order to characterize in probabilistic terms the mean velocity of a probable wind climate (Torrielli et al. 2010). Once such a stochastic process is fully described in terms of its marginal probability distribution function (PDF) and power spectral density function (PSDF), it is simulated by means of a suitable algorithm, recently proposed by the authors (Torrielli et al., 2010; 2011), that is a modified version of the Masters and Gurley’s (2003) method based on the theory of translation processes (Grigoriu, 1995). Through the application of this type of procedure, 386 33-year time series of mean wind velocity are generated, forming a set of over 12,000 years of synthetic wind speeds. EVs from this vast dataset are fitted to the EV distributions discussed in Section 2, so as to assess their suitability for use in modeling the annual maximum speed distribution. This large quantity of data reduces the negative effect of sampling errors, thus allowing a more accurate evaluation of the distributions investigated. In this regard, a parametric analysis is also performed to study how record length affects the effectiveness of the EV models and to investigate their applicability to real cases. Section 5 provides overall considerations about the EV distributions used to model the traditional dataset of wind speeds, taking into account the theoretical issues, the findings resulting from the analysis of the large sized dataset, and the applicability of each EV model to real or probable situations (Section 4).

2. Extreme wind distributions The following sub-sections provide a synthetic framework of the extreme wind distributions investigated and compared in the present paper. 2.1. Population analysis Classical probability theory states that the largest of n independent samples, V^ n, extracted from a set of data of the same parent distribution FV(v), is distributed according to the first of Eq. (1). PrfV^ n o vg ¼ ½F V ðvÞn ;

PrfV^ T ovg ¼ ½F V ðvÞrT ;

rT rn

ð1Þ

When a serially-correlated stochastic process is considered, the hypothesis of independent samples is no longer applicable. In this case, the distribution of V^ T , the largest value in an epoch of T years, can still be expressed by the preceding formula, provided it is corrected as shown in the second of Eq. (1), where the exponent n is replaced by rT (Davenport, 1968; Cook, 1982), r being the annual rate of independent events. Unfortunately, commonly-applied methods to estimate r are rather unreliable, for this reason asymptotic analysis being preferred since it assumes rT-N.

23

2.2. Asymptotic analysis Focusing on the second of Eq. (1), Fisher and Tippett (1928) showed that if rT approaches infinity, then the cumulative distribution function (CDF) of the largest value in period T tends to one of three possible distributions. Alternative derivations of the main results of EV analysis, provided by Cramer (1946), give a particularly clear understanding of the basics of asymptotic analysis. In this regard, the second of Eq. (1) can be written:  rT   Q V ðvÞ Q ðvÞ PrfV^ T ovg ¼ ½F V ðvÞrT ¼ 1 ð2Þ !exp  V rTUQ V ðuT Þ rT-1 Q V ðuT Þ where QV ¼1 FV is commonly called risk of exceedance, while the characteristic largest value uT for a period T is the value that meets the condition QV(uT) ¼(rT)  1. For all v4uT, QV(v)/QV(uT)o1 and consequently if rT is large enough, assuming rT-N, the Cauchy (1821) formula provides the RHS of Eq. (2). It is worth noting that this is a universal result regardless of the shape of the parent distribution. Hereafter, modeling of the ratio QV(v)/QV(uT) depends on the shape of the parent and determines the asymptotic form. The three asymptotes, named later by Gumbel (1958) as Fisher and Tippett Type I (FT1), Type II (FT2) and Type III (FT3), had been previously combined into a single distribution by Von Mises (1936, in French; for English see Jenkinson, 1955), now known as the generalized extreme value (GEV) distribution: (   ½1b aðvuT Þ1=b b a 0 PrfV^ T ovg ¼ exp LðvÞ ; LðvÞ ¼ ð3Þ exp½aðvuT Þb ¼ 0 where uT and a  1 are, respectively, the location (mode) and dispersion (or scale) parameters, while b is the shape factor which determines the asymptotic form assumed by the GEV distribution; b o0 corresponds to the FT2 or Frechet distribution, b ¼ 0 corresponds to the FT1 or Gumbel distribution, and b 40 corresponds to the FT3 or reverse Weibull distribution. Several techniques to estimate the GEV parameters exist, the most common described in literature being probability weighted moments (PWM) and maximum likelihood (ML) (e.g. see Palutikof et al., 1999). However, as pointed out by Harris (2006), no data fitting could yield the condition b ¼0, since it is associated with a singularity of the exponent. This fact renders necessary the use of an effective test to verify the shape factor estimated. For this purpose, Hosking et al. (1984) furnishes a simple but effective goodness-of-fit test called the Z-test: rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi N Z¼b ð4Þ 0:5633 where the Z-statistic is normally distributed with zero mean and unit variance, while N is the size of the sample analyzed to evaluate b. Z is compared with reference values associated with prefixed confidence limits to establish if the b estimation is correct or an artifact of sampling error. 2.3. FT1 penultimate distribution for tail-equivalent Weibull parents Asymptotic distributions are often applied in spite of the fact that rT is too low to meet the basic assumption rT-N. In fact, in wind engineering design, typically T¼1 year and practical values of r may be 876 (estimated at Brookhaven, USA, by Davenport, 1968) or 797 (estimated at Boscombe Down, UK, by Cook and Harris, 2008). Consequently, it is common that the GEV distribution results in an asymptote incompatible with the upper tail of the parent (Harris, 2004). To avoid this situation, Cook and Harris (2004, 2008) proposed the use of FT1 penultimate distribution, which takes into account the actual finite values of rT. It assumes

24

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

a very elegant shape when extremes are drawn from a tailequivalent Weibull parent: n h

io ; PrfV^ T o vg ¼ exp exp ak vk ukT

a¼

1 c

ð5Þ

2.4. Process or level-crossing analysis The idea underlying this method was proposed by Davenport (1968) in order to overcome the hypothesis of independent events forming the basis of EV methods. Gomes and Vickery (1977) restored the basic idea of considering wind speed as a stationary random function; according to this hypothesis, the annual mean number of up-crossings of a threshold v can be expressed by (Rice, 1944, 1945): Z 1 þ _ V V_ ðv, v_ Þdv_ NV ðvÞ ¼ ð6Þ vf 0

where V_ is the derivative of the process V, and f V V_ ðv, v_ Þ is the joint PDF of V and V_ . When the velocity thresholds are sufficiently high, the related up-crossings can be regarded as rare and independent events, and thus treated as a Poisson process. Under such an assumption, the CDF of V^ T can be identified with the probability of zero crossings of level v in epoch T, and can be expressed by: þ

ð10Þ þ

where c and k are the parameters of the Weibull parent, while a and uT are related to the EV distribution. If the population distribution is known, only uT, or alternatively both a and uT, are estimated from the extremes; when the parent is unknown, k becomes an additional degree of freedom and a three-parameter fit is performed by means of an iterative procedure (Cook and Harris, 2004). Generally-speaking, it is reasonable to expect the fit to the population to lead to more reliable estimations than the fit to the annual maxima, as the larger number of observations narrows the confidence limits of the estimated parameters. From this point of view, the FT1 penultimate estimating only uT mode is to be preferred when the parent is known; however use of the other applications of the FT1 penultimate is justified when the parent parameters are not estimated accurately.

PrfV^ T o vg ¼ exp½TUNV ðvÞ

Vickery (1977): " þ # N ðvÞ l¼E V f V ðvÞ

ð7Þ

According to the multiplication theorem of probability, f V V_ ðv, v_ Þ can be expressed as the product of the marginal PDF of V, fV (v), and the conditional PDF of V_ given V, f V V_ ðv, v_ Þ. Hence, alternatively, Eq. (6) can be rewritten as follows: Z 1  þ _ V_ 9V v9v _ vf dv_ ð8Þ NV ðvÞ ¼ lðvÞf V ðvÞ; lðvÞ ¼ 0

When V and V_ are statistically independent, f V V_ ðv, v_ Þ can be expressed as the product of the marginal PDFs of V and V_ . Although this is true if V has a Gaussian or Rayleigh-distribution, this assumption was considered plausible by many researchers, including Davenport (1968). Under such a condition, Eq. (8) becomes: Z 1 þ _ V_ ðv_ Þdv_ vf ð9Þ NV ðvÞ ¼ l f V ðvÞ; l ¼ 0

where the overhead symbol is used to emphasize that l is independent of v. Even if Eqs. (8) and (9) show strict formal analogies, it is worth noting that the former rejects, while the latter accepts, statistical independence between V and V_ . The evaluation of l as defined in Eq. (9) is not straightforward since it requires knowledge of f V_ ðv_ Þ that is normally burdensome to evaluate. An easier technique was suggested by Gomes and

where E[  ] is the mean operator; N V (v) is numerically estimated by enumerating all up-crossings of level v in epoch T; and fV (v) is the PDF of the parent population. Unlike previous methods, the process analysis does not describe the EV distribution by directly analyzing the extremes, but by analyzing the parent. Therefore, even short wind records allow the related EV distribution to be described. 2.5. Peak over threshold (POT) method The POT approach is the only method that allows for both the use of sub-annual maxima and the option of a GEV fit. Pickands (1975) considered a variable V having a parent FV(v), such that, for large rT, the distribution of the largest value in the period T converges to one of the asymptotes combined in Eq. (3). He showed that the excesses of such variable over a threshold u, X ¼ Vu, tend to conform to the Generalized Pareto Distribution (GPD): h x i1=b F X ðxÞ ¼ 1 1b ba0 d ð11Þ x

b¼0 F X ðxÞ ¼ 1exp 

d

provided that u is arbitrarily large. In Eq. (10) d is the scaling parameter and b is the shape factor previously introduced for GEV analysis. Moreover, if the threshold is large enough, its crossings can be assumed to be independent, and the number N of values over the threshold u in a period T is Poisson-distributed, with a rate parameter lu/year. An unbiased estimate of lu is n/T (Abild et al., 1992), where n is the total number of exceedances of u counted directly from the data, and T is the number of years of the record. In such a case, the distribution of the largest value of V in T, V^ T, is given by (Abild et al., 1992): (   ) ðvuÞ 1=b ^ PrfV T o vg ¼ exp lu T 1b b a0

d

  ðvuÞ PrfV^ T ovg ¼ exp lu T exp 

d

b¼0

ð12Þ

From Eq. (12) it appears that the distribution of V^ T is highly dependent on the threshold u, both directly and indirectly, through lu, d and b. On the one hand, u has to be set low enough to ensure that a sufficient quantity of data is used to estimate the distribution parameters, while on the other hand the asymptotic requirement u-N, supporting the GPD distribution, has to be satisfied. Furthermore, a suitable combination of the threshold u and the minimum separation time between events is required, so as to guarantee that only independent exceedances are counted in the estimation of lu (Walshaw, 1994). Since u is chosen by the user and lu is derived from the data, the POT method reduces the fitting problem to estimating only d and b instead of the three parameters required by the GEV approach. Many methods are available to fit data sets to the GPD model in Eq. (11), such as maximum likelihood (ML; Davison and Smith, 1990), probability weighted moments (PWM; Harris, 2005) and conditional mean exceedances (CME; Davison 1984). 2.6. r-Los method This approach selects the r largest order statistics (r-LOS) in a single epoch, usually one year, to estimate the parameters of the GEV distribution. The theory behind this method was developed

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

by Weissman (1978), but applications in wind engineering are quite rare. The method is based on a joint distribution of a single r-LOS vector (Xr oXr-1,y,oX1) that can be derived using the theory of Poisson processes (Coles, 2001) as: f 1,::,r ðx1 ,::,xr Þ ¼ exp½Lðxr Þ

r Y

a½Lðxk Þ1b

ð13Þ

k¼1

where a, b and L(x) are defined in Eq. (3). This joint distribution provides the basis for the ML method. Suppose that R years of data are available and r-LOS values are extracted from each year of data, the likelihood function can be obtained as a product of the R densities given in Eq. (13) corresponding to the observed sample: ( ) R r Y Y LðuT ,a, bÞ ¼ exp½Lðxr,i Þ a½Lðxk,i Þ1b ð14Þ i¼1

k¼1

where xk,i denotes the k-th largest OS in i-th year of the record. The parameters of the GEV distribution are finally derived by maximizing the log-likelihood function. The selected r-LOS values must be independent events, in line with the basic assumption for the EV analysis. It is reasonable to assume that wind speeds from different years meet this criteria; on the other hand, the strong serial correlation of the velocity process requires that a suitable separation interval be set between the data of the same year. This issue is strictly related to the choice of the number of LOS involved in the analysis, r (Smith, 1986). A practical criterion is to set r so that it minimizes the variance associated with the quantity to be estimated. In this regard, in previous applications of the r-LOS method the setting was r ¼3 to r ¼7 by Tawn (1988), r ¼10 by Coles and Walshaw (1994), and r ¼5 by An and Pandey (2007). In conclusion, it is reasonable to regard the r-LOS method as an extension of the classical asymptotic analysis, where more values are extracted from each year of observation. In this regard, Harris (2009) made clear that this method also rests on the asymptotic assumption of rT-N, which is not satisfied considering the actual number of independent events in T; therefore the r-LOS method is subject to the same systematic bias errors as other GEV methods. 2.7. Methods based on independent storms (MIS, IMIS and XIMIS) The Method of Independent Storms (MIS) was originally devised by Cook (1982), and subsequently improved by Harris (1999) becoming known as IMIS. It examines continuous wind records in order to detect independent storms and their maximum values, with a typical rate rs of approximately 100 events/year. Firstly, a filter is applied for the purpose of removing peaks of short-duration caused by meteorological mechanisms other than the dominant mechanism (Cook, 1982); it consists in the calculation of 10 h non-overlapping moving averages. Secondly, the filtered time series is searched to find the lulls, identified by downward crossings of a prefixed low threshold. Between each pair of lulls there is a storm, whose independence from the preceding and succeeding events is guaranteed by the presence of the lull; the maximum speed occurring during each storm is established considering the original unfiltered record. After detecting the extreme velocities, they are transformed into dynamic pressures to achieve a more rapid rate of convergence, and the top M order statistics are fitted to the Type I asymptote. Harris (1999) highlighted the systematic errors in Cook’s original method induced by both the estimator employed for the annual maximum speed, and the method used for fitting the Type I asymptote. To eliminate these errors, Harris proposed a new set of plotting positions based on the mean values of the order statistics.

25

He further replaced the original Lieblein BLUE fitting technique (Lieblein, 1974; Cook, 1985) with a new minimum-variance biasfree procedure (Harris, 1996), based on the weighted least squares (WLS), involving the standard deviations of the order statistics in the weight definitions. Unfortunately, as Harris (1999) himself recognized, the factorial function prevents the applicability of this procedure in the case of samples containing a large number of extremes. This limitation is removed by the XIMIS method, introduced by Harris (2009) as a further refinement of IMIS, and suitable for all types of wind climate. Some improvements have been made to the modeling of the distribution of the largest values over the period T; in particular, the standard reduced variate of V^ T : y ¼ ln½lnðPrfV^ T o vgÞ

ð15Þ

is modeled by a new plotting position: ym ¼ lnðRÞcðmÞ

ð16Þ

where R is the number of complete epochs T covered by the speed record; c(v) is the Psi Function, i.e. the logarithmic derivative of the Gamma Function; and m is the rank. The plotting positions are assigned to the data sample ranked in descending order, so m¼1 is associated with the largest value. It must be stressed that the formulae of the plotting positions in Eq. (16) can be applied only if the dataset is large enough, say N Z30, and it has the advantage of not requiring the annual rate of independent events, r, which is hard to define. Furthermore, since the basic assumption of asymptotic analysis rT-N is not made, this method is ascribable to the penultimate distribution family. A distinguishing feature of the XIMIS method is that prediction of standard design values, such as the 1:50 yr value, does not require any sort of probability model, but relies solely on experimental evidence. In particular, practical applications suggest the existence of a power law transformation (Naess, 1998; Harris, 1999; Palutikof et al., 1999) Z ¼Vw such that the resulting data values, ranked in descending order Zm and plotted against the reduced variate ym , are well approximated by a straight line: ym ¼ aZ m P

ð17Þ

The power w has to be chosen to produce the best fit of the data to the straight line, with minimum residual. The slope and intercept of the resulting straight line define, respectively, the a and P parameters in Eq. (17), knowledge of which allows V50 to be estimated through interpolation or, when appropriate, mild extrapolation. The prediction of extreme values with a high return period, often assumed as the design values for strategic facilities, cannot rely on simple empirical considerations, as their values lie well beyond the range of commonly-available data. In this case, it is necessary to define in advance an EV model, Harris suggesting the use of the Type I asymptote as the basis for extrapolation. It is worth noting that, when the parent is Weibull or a distribution of exponential type, a linear dependence between Zm and ym is assumed, and the prediction is again reduced to a simple linear extrapolation, similar to the previous estimate. In the light of the previous discussion about XIMIS, considering the transformation Z¼Vw and Eqs. (15) and (17), it is trivial to demonstrate that the distribution of the annual maximum velocity can be defined by the following expression:  PrfV^ T ovg ¼ expfexp½ avw P g ð18Þ This type of distribution is independent of both the number of independent events rT and the parent distribution. It is interesting to observe the close relationship between XIMIS and FT1 penultimate distribution in Eq. (5). Even in this case rT is not required, and knowledge of the parent is not strictly necessary as long as it

26

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

is known to be tail-equivalent Weibull. In addition, both distributions have the same double exponential form, and employ preconditioning to obtain a better convergence with the asymptotic form. In particular, when the parent distribution is not available, the estimate of the FT1 penultimate parameters k, a and u is totally linked to the extremes, the same being true for the XIMIS parameters w, a and P. In this case, if the fitting procedure used is the same, the XIMIS and FT1 penultimate models resemble each other. Cook and Harris (2004) for the FT1 penultimate, and Harris (2009) for the XIMIS, proposed two different techniques to perform a three-parameter fitting, both aimed at minimizing the residual of the WLS method proposed.

3. Simulation of mean wind speed process The wind speed process V(t) to be simulated is described in terms of a probability distribution function (CDF) and a power spectral density function (PSDF) derived from real meteorological data. 3.1. Wind speed databases This paper resumes the case study extensively presented in Torrielli et al. (2010). The description of V(t) is based on the wind measurements recorded by 6 anemometers spread out over an area of approximately 200 km in diameter in the central part of Italy. Such measurements comprehend both historical and modern records. The historical records, provided by the Air Force and Ente Nazionale per l’Assistenza al Volo (Flight Control Agency), are a collection of speeds taken over a 33 year period at 3 h intervals. The modern records, provided by the new instrumental network of Rete Ferroviaria Italian (Italian Railway Company), are a collection of speeds taken over a 1 year period at 10 min intervals. In both cases, the velocity values are averaged over 10-min intervals. Three sources of possible inconsistencies – namely wind calms, local terrain and mixed-climate features – are firstly analyzed and mitigated as much as possible. The problem of wind calms results from the replacing of obsolete instruments and the increasing sensitivity of modern ones. This was highlighted by the fact that recent wind calms are much less frequent than in the past. Thus, older wind calms clearly include false observations that may prejudice the representativeness of the older records. To limit this problem, and render wind records more homogeneous, the most recent and reliable measurements are used to correct probabilistically older, less accurate ones (Torrielli, 2011). The features of the local terrain where the 6 anemometers are sited are taken into account by transforming the wind data – affected by the roughness and topography of the local terrain – into speed values corresponding to a unique reference condition, i.e. at a height of 10 m on a flat terrain with a roughness length z0 ¼0.05 m (Burlando et al., 2010). The fulfillment of this requirement is made considerably simpler because all sensors were placed in flat open territories, sufficiently far from mountains and the sea coast. This also reduces the possibility that local up/downhill winds or breezes compromise greatly data homogeneity. Mixed-climate features that could cause extreme wind speeds associated with different meteorological phenomena – e.g. extratropical cyclones on a synoptic scale and thunderstorm winds on a mesoscale – represent the most delicate issue. In principle, wind speeds originating from different types of phenomena should be separated carefully and dealt with independently. However, this type of operation is extremely difficult and with scarce reliability in all cases where only average and peak values are available at 10 min/1 h intervals (in practice excluding all case studies). In

addition this operation is almost impossible in Italy, where historical records are related to wind speeds 3 h apart, averaged over 10-min intervals, and peak wind speeds have been recorded only during the most recent years. In spite of this, a significant effort has been made to remove isolated, intense wind speeds not ascribable to synoptic scale cyclonic events from the database. Thus all databases are rendered as homogeneous as possible and representative, in terms of extremes, of a unique wind type. A modern anemometer network has been recently created by two of the authors in the High Tyrrhenian area (Solari et al., 2012), consisting of 33 ultra-sonic anemometers and an acquisition system for the purpose of creating a huge database of wind speeds continuously recorded as time steps of 0.1 s. A wide research program is underway to establish an efficient procedure for identifying and systematically separating wind speeds associated with different events. 3.2. Mean wind speed process The 6 wind speed records described in the previous section have been processed in order to characterize a wind speed process V(t) associated with a probable wind climate. This definition means that although the latter is not retraceable to a given geographical site, it is nevertheless physically plausible for the area investigated. In particular, the CDF of V(t) is modeled by considering only the 3 historical records; collecting 33 year of measurements, these are more representative than the 3 recent records. Conversely, the definition of the PSDF uses both historical and modern speed records. The Weibull model is the most commonly used to represent the parent distribution of the mean wind speed; however, it does not consider the presence of wind calms. Even though the actual existence of wind calms is a controversial matter, surely the time correlation of the current V(t) has to be computed considering a record with zero values. So, the presence of wind calms should be taken into account also in the probability distribution of the parent, at least to be consistent with the correlation of the process. In light of these considerations, the marginal CDF of V(t) is described by the Hybrid Weibull (HW) model (Takle and Brown, 1978), which corrects the Weibull distribution to consider the presence of wind calms:

 k  v F V ðvÞ ¼ F 0 þ ð1F 0 Þ 1exp  ð19Þ c where F0 is the rate of wind calms, while c and k are the Weibull parameters. It is trivial to demonstrate that the model in Eq. (19) is a tail-equivalent Weibull distribution. The original 3 datasets of historical wind speeds, suitably homogenized, are merged into a single database associated with a so-called superstation (Peterka, 1992; Simiu and Filliben, 1999). Such a dataset is fitted to the HW parent distribution so as to estimate the related parameters F0, c and k reported in Table 1 together with the first four central moments, i.e. mean mV, variance sV 2, skewness g1 and kurtosis g2. Note that the superstation dataset is studied only to describe the parent distribution, and the EV distribution is not considered. It is worth noting that the value assumed by the shape parameter k¼1.16 is quite outside the range 1.8–2 found in UK or in other Northern European countries (Cook and Harris, 2004). Table 1 Parameters of HW distribution and the relative moments. F0

k

c (m s  1)

mV (m s  1)

s2V (m  2 s  2)

g1

g2

0.118

1.155

3.091

2.592

6.637

1.611

6.631

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

Nonetheless, it is well in line with previous estimates related to the wind climate prevalent in the central regions of Italy (Freda and Solari, 2010). The decision to develop the present example for a wind speed process with so small a k value is based on using the most reliable data from those currently available. Much improved data is currently being collected in Ligury (Solari et al., 2012), where k¼1.4–1.7. As soon as the quantity of data gathered is large enough, a new wind speed process will be examined to obtain more general results from this research, and to discuss the role of k. The possibility of generating numerically a catalog of wind speed processes coherent with target analytical PDFs and PSDFs representative of a broad band of worldwide wind climate conditions is also being considered. Taken separately, the historical and modern records are not able to describe the macro-meteorological spectral component for the whole frequency range of interest. The combined use of these records ensures a complete spectral description of the mean wind speed. This operation is performed by the spectral combination technique presented in Torrielli et al. (2010, 2011). Firstly, the classical spectral analysis is employed to estimate the PSDFs of the 6 anemometer records. The spectra of the historical records cover periods from 6 h to 33 years, while those of the modern records cover time intervals from 20 min to 1 year. Next, the mean spectra associated with the historical (S1) and modern (S2) records (Fig. 1b and c) are derived by averaging the PSDFs of each group. The derivation of the total spectrum of the mean wind speed is straightforward for period/frequency bands where only a source spectrum is defined, i.e. S1 or S2. In the overlap region, the total spectrum is obtained by averaging S1 and S2, weighted, respectively, by the cosinusoidal functions W1 and W2 (Fig. 1a and d), equal to half Hanning windows. Their shapes are particularly advantageous: W1 reduces the contribution of S1 for periods of approximately 6 h, where spectral estimates are corrupted by the aliasing effect; W2 fades the influence of S2 for approximately 1 year, where estimates are poor because records contain few complete oscillations of such spectral components. Fig. 1e shows the resulting PSDF of the macro-meteorological component of the velocity of the probable wind climate, suitably normalized to remove the alterations of energy/variance induced by the weighted average. It presents a series of spikes corresponding to the annual cycle (1 year), the diurnal cycle (24 h) and its super-harmonics (12 h and 8 h). The decreasing trend shown on

27

the right side of Fig. 1e reflects the spectral gap. These features are very similar to other horizontal surface wind spectra reported in literature, such as the spectrum at Matsumoto station analyzed by Kai (1987). It is worth noting that Harris (2008) points out that the annual cycle, diurnal cycle and its super-harmonics are periodic components, with a fixed phase associated with the local wind climate. These periodic components are usually subsumed into the mean wind speed record when extreme value analyses are carried out. This might increase the dispersion parameter of the extreme value distribution, since from 1 year to another the true maximum will not coincide with the same portion of any periodic cycle. In the present paper such components are preliminarily involved in the simulation, and treated as stochastic terms with random phase. A new paper to be submitted for publication by these same authors (Torrielli et al., in preparation) analyses the effect of the deterministic components, demonstrating that the approximation in the present paper does not compromise the independence of the extreme values of the simulated timehistories. 3.3. Wind velocity simulation The mean wind velocity process is generated by means of a non-Gaussian simulation method based on memoryless transformations, which has recently been proposed by Torrielli et al. (2011). This algorithm corrects iteratively the probabilistic and spectral content of the sample function to be generated by applying, respectively, a CDF map-based transformation and the spectral correction technique (Fig. 2). At the end of each iteration, the spectral content of the sample function v(t) matches the target spectrum, but its distribution is distorted from the target distribution. The probabilistic distortion induced on v(t) is evaluated by comparing the sample skewness and kurtosis with the corresponding target values. If the differences are below prefixed tolerance values, the algorithm stops, otherwise a further iteration is performed. The algorithm, herein applied for simulating HW-distributed processes, is a revised version of Masters and Gurley’s (2003) method. Masters and Gurley’s method corrects iteratively the probability distribution of a given sample function v(t) by means of memoryless transformations. To achieve this, a semi-empirical model is used to mimic the unknown distribution of v(t).

frequency [Hz] 1

6

W1

0.5 0

-6

10

-4

24 hr

5

nS (n) [m2/s2]

5

S2 10

10

5

S1 10

10

-8

4 3 12 hr

2 8 hr

1

1 W2

0.5 0 10

33 yr

4

1 yr

10

period [hr]

2

10

6 hr

0

20 min

1 yr

0 10

4

10

2

10

0

period [hr]

Fig. 1. (a–d) Sketch of the procedure used to set the PSDF of the probable wind climate; and (e) normalized PSDF of the mean speed of the probable wind climate.

28

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

Fig. 2. Algorithm used for the simulation of the mean wind speed process.

7 simulation target

autocorrelation (m2/s2)

6 5 4 3 2 1 0 -1 0

20

40

60 80 time lag (hours)

100

120

140

Fig. 3. Autocorrelation of the simulated mean wind speed process.

The original Masters and Gurley’s method guarantees that only the maximum and minimum values are random. The goal of the modified algorithm is to generate sample functions having more extreme values that are all really random. For this purpose, this paper replaces the empirical model introduced by Masters and Gurley with a semi-empirical model. The semi-empirical model consists again of an empirical model (based, for instance, on given plotting positions) for describing the unknown distribution of v(t) for negative values, while a HW distribution is used for non-negative values. The parameters of the HW distribution have to be estimated for each iteration. With reference to the j-th iteration, the distribution of sample vj(t) is given by:  8 vj o 0 > F V j ðvÞ  <

  

½kj  F V j vj ¼ ð20Þ vj > vj Z 0 : ½F 0 j þ 1½F 0 j 1exp  ½c j

where F V j ðvÞ represents the empirical CDF of vj(t), and [F0]j, [c]j, [k]j are the HW parameters associated with the j-th iteration. The use of an analytical model for describing the positive values, and hence the upper tail of the distribution of v(t), allows sample functions to be generated with highest values independent of sample size and so effectively random. For a more detailed description of the simulation algorithm the reader is referred to Torrielli et al. (2011). This simulation algorithm is applied to generate 386 time series of the mean velocity of the probable wind climate. Each time series consists of 33 years of 10 min average wind speeds at 10 min intervals. Note that the simulated wind speeds are in m/s and are not rounded to the nearest knot. This is done to avoid rounding errors which usually affect real wind data.

generation of sample functions that match precisely the prescribed PSDF. The related marginal CDF is slightly distorted from the target distribution, coherent with the adopted simulation method. It is worth noting that the empirical distribution of the synthetic sample function is modeled here by the order statistic median, instead of the more commonly-applied order statistic mean. Although analytical formulae for the plotting position of the order statistic median of N correlated observations are available, for instance Young (1967), their evaluation is very burdensome for the current value of N ¼1735761. As an alternative, the plotting positions of the order statistic median can be estimated through the set of the simulated sample functions, consisting of 386 items. In particular, the distribution of a generic wind speed function vr(t) is modeled by the pairs of points (vr(i), Med{F[vr(i)]}), where {vr(i):vr(1) ovr(2) oyovr(N)} is the sample consisting of the vr(t) values ranked in ascending order, F is the target CDF, and Med{  } is the median operator applied to the i-th set {F[vr(i)]} with r ¼1,...,386. These plotting positions are used in Fig. 4a to model the distribution FV of a single simulated sample function ( þ), while the straight line (–) is the target CDF. The probability paper shown in Fig. 4 is a modified version of the classical Weibull plot to represent the current Hybrid Weibull parent. This diagram provides ln(v) on the abscissa and ln{  ln[(1 FV)/(1 F0)]} on the ordinate, where ln{  } is the natural logarithm of {  }. These transformations of the axes were chosen to ensure that the CDF of a HW random variable results in a straight line if plotted on probability paper. Fig. 4a shows good agreement between the sample and target CDF over the whole range of velocities, with the exception of the upper tail where a slight detachment occurs. Two factors contribute to this discrepancy: the statistical uncertainty deriving from the generation of finite-size (N) sample functions, and the proposed simulation algorithm which is not able to reproduce exactly the desired distribution. To investigate the nature of these distortions, Fig. 4b does not show the single but the ensemble distribution relative to the 386 sample functions (þ). In this case, the empirical CDF is modeled through the pairs of points (Med{vr(i)}, Med{F[vr(i)]}). The perfect agreement between the empirical and target distributions reveals that the distortions associated with the various sample functions offset each other, as in ensemble statistics. This observation implies that the numerical regression of the EV samples induces a mutual compensation of those distortions, and supports analysis of the extreme values deriving from the simulated sample functions. Various studies are currently being carried out by the authors to overcome the issue of probabilistic distortions by reproducing the same analyses with the use of different simulation algorithms.

3.4. Examination of simulated sample functions 4. Extreme value analysis of long-term wind time series Fig. 3 shows the autocorrelation functions of the target process V(t) (—) and of a generic simulated time series (þ). This demonstrates that the algorithm applied guarantees

The simulation of the mean wind speed process in Section 3 provides 386 33 year-long time series. Thus, a large sample of

ln(-ln[(1-FV)/(1-F0)])

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

1-1E-7

1-1E-7

1-1E-5 1-1E-3

1-1E-5 1-1E-3

1-1E-2

1-1E-2

0.9

0.9

0.8 0.7 0.6 0.5 0.4

0.8 0.7 0.6 0.5 0.4 2

5

10

15 20 253035

ln(v)

29

Target emp CDF 2

5

10

15 20 253035

ln(v)

Fig. 4. Plotting positions of the order statistic median: (a) single sample functions; and (b) ensemble distribution.

Fig. 5. Gumbel plot for annual maximum mean wind speeds from annual maxima.

wind speeds consisting of 12,738 annual maxima is investigated by applying the EV techniques introduced in Section 2. 4.1. Annual maxima The first step towards a direct investigation of the annual maximum values is to detect the largest speeds at 1-year time intervals, from each of the simulated time series; these values are then gathered into a single large sized dataset consisting of M¼12,738 annual maxima. As pointed out by Harris (2004, 2009), if rT is large enough (say Z50), as is usual for temperate storms, all EV methods based on the Cauchy formula (such as FT1 ultimate, GEV, FT1 Penultimate, POT, and XIMIS) can be applied provided that the annual maxima are greater than a characteristic largest value u1. Harris (2004) outlines a simple method for estimating the characteristic largest value, independently of which asymptote is fitted. In the present study u1 ¼22 m/s is assumed to be the lowest allowable speed value, and annual maxima below this limit are rejected for the analysis. In this way, the original dataset consisting of 12,738 items is left-censored, and replaced by a new smaller dataset consisting of 7991 annual maxima. The related sample distribution () is plotted in Fig. 5 on Gumbel probability paper. The bootstrap procedure described in

Cook (2004) is applied to derive the mean plotting positions for representing the empirical distribution, and the weights for the WLS method proposed by Harris (1996). Alternatively, the plotting position formulae of the order statistic median may be applied, as those proposed by Yu and Huang (2001) for uncorrelated observations. The empirical distribution () shows a regular tendency, with upward concavity, in the range of mean wind velocity vr35 m/s. Above this level some fluctuations occur, and the empirical distribution loses its reliability due to the finite size of the dataset analyzed. The large sized dataset of annual maximum speeds is fitted to three EV distributions: the asymptotic FT1 ultimate and GEV distributions in Fig. 5a, and the FT1 penultimate distribution in Fig. 5b. Three different sub-cases for the FT1 penultimate distribution are investigated, in accordance with Cook and Harris (2004): case (1), the parent is unknown and so k, a and u1 must be evaluated by means of a three-parameter fitting; case (2), the shape parameter k of the parent is fixed, while the estimation of a and u1 is regressed to the extremes; and case (3), the parent is known and so only u1 needs to be estimated from the data. Harris’ WLS method is applied to estimate the parameters of FT1 ultimate, penultimate and GEV distributions. Note that the PWMs method, applied commonly to estimate the parameters of

30

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

Table 2 Parameters of the annual maximum speed distributions from the annual maxima. Values in bold do not result from the regression of the extremes but are assigned in advance. FT1 ultimate

GEV

b a [m  1 s] u1 [m s  1]

0.404 22.04

r

17,953

a u1 L [m s  1] r

0.074 0.369 21.96 58.58 17,231

FT1 penultimate (fit k, a and u1)

FT1 penultimate (fit a and u1)

FT1 penultimate (fit only u1)

k a u1

1.680 0.117 21.98

k a u1

1.155 0.262 22.03

k a u1

1.155 0.324 22.41

r

130

r

15,726

r

19,070

the GEV distribution, does not support left-censoring. For this reason, use of the WLS method is preferred. Some directions on the use of the WLS method for estimating the parameters of the GEV distribution are given by Cook and Harris (2004). Table 2 lists the resulting estimates of the parameters. It is worth noting that the mode u1 of annual maxima is very consistent across EV models. Bearing in mind the definition of uT in Eq. (2), comparable estimates of the annual rate of independent events r can be derived. The only exception is given by the three-parameter fit of the FT1 penultimate. On the other hand, the estimates of the dispersion parameter 1/a cannot be directly compared, since its definition changes from method to method. Similar results were obtained by Cook and Harris (2004) analyzing 30 years of wind observations at Boscombe Down. Fig. 5a illustrates clearly that the straight line of FT1 ultimate distribution fails to capture the upward concavity of the annual maxima at the upper tail of the empirical distribution, for say, v430 m/s. Hence the FT1 ultimate distribution yields to quite conservative design speeds if they are associated with high return periods. The GEV distribution results in the FT3 asymptotic form, since b ¼0.062 (Table 2) is greater than 0. This finding is verified by applying the goodness-of-fit Z-test proposed by Hosking et al. (1984). All significant tests, including the Z-test, deal only with random sampling errors, presuming that all systematic errors have been eliminated. A preliminary estimate of the systematic error e 0.02 can be derived from the error chart provided by Harris (2006), considering rE 17,231 (from Table 2) and k¼1.16, respectively, the rate of the independent events and the shape parameter of the HW parent known in advance. The systematic error is removed from b and the resulting value is applied to evaluate the Z statistics by means of Eq. (4). The result Z ¼5.00 is greater than the value required for the significance test at the 5% level, i.e. 1.645, and confirms the FT3 asymptotic form. Although the FT3 asymptotic form is confirmed by the significance test, the incompatibility between the FT3 asymptote and the exponential parent from which the maxima are drawn should not be overlooked. On the other hand, Fig. 5a shows that the GEV/FT3 distribution fits the data very well, at least for vo40 m/s. It is difficult to judge the fitting of the upper tail due to the fluctuations of the sampling data. The single parameter FT1 penultimate (y) distribution in Fig. 5b tends to underestimate the design speeds associated with high return periods. Conversely, the FT1 penultimate two-parameter distribution (- -) seems to overestimate the design speeds associated with high return periods. The three-parameter FT1 penultimate (–) yields to a careful modeling of the distribution of the annual maxima, in the whole range of wind speeds. Unfortunately, the k and a parameters (Table 2) derived by fitting the extreme values to the three-parameter FT1 penultimate are not consistent with the parent parameters (Table 1), which were assigned to the time series during the simulation. This highlights a lack of agreement between the true parent and the parent distribution description derived by processing the

Table 3 Parameters of the process analysis distributions. Grazzanise Pr1a (v40)

Probable wind climate Pr1a (v40) Pr1b (v425)

Pr2

l [y  1]

l

l

d [yr  1  m  1  s]

g

2853

41145

70592

181.48

1.762

extreme values drawn from the true parent. Cook and Harris (2004) found the same incompatibility analyzing 30 years’ wind observations at Boscombe Down. They identified the causes of such disagreement in the description of the population by means of a single-mechanism Weibull model (Cook and Harris, 2008). However, it seems unlikely that the same causes may explain the inconsistency that occurs in this case, where a single-mechanism climate is expected (Section 3.1). 4.2. Process or level-crossing analysis The process or level-crossing analysis does not involve any asymptotic assumptions (Cauchy formula) so, generally, data do not require a left-censor. In order to ensure a correct application of this method, two key issues have to be considered: (1) the assumption of independent up-crossings which supports the þ Poisson process; and (2) modeling of NV (v) in Eq. (7), i.e. the annual mean number of up-crossings of threshold v. Both these þ issues can be addressed by estimating NV (v) numerically from the simulated wind speed time histories. The opportunity to analyze more than 12,000 years of wind speeds is unprecedented, þ and allows a very accurate description of NV (v). Three distinct applications of the process analysis, namely Pr0, Pr1 and Pr2, are presented in this study. þ Pr0 is based on the numerical estimate of NV (v), directly determined by counting the up-crossings of the simulated time histories. On the one hand, this makes Pr0 a rigorous application of the process analysis while, on the other, it cannot be used to derive design speeds associated with return periods longer than the record length. Conversely, other applications of process analysis involve analytical methods, each characterized by different modeling of þ NV (v). Pr1 is the classical application of process analysis, based on þ the assumption of independence between V and V_ . N V (v) is estimated by means of Eq. (9), where the parent distribution is known in advance in this study; l is usually estimated by Eq. (10), þ where N V (v) is numerically determined, and fV(v) is the PDF of the parent distribution. Two cases are considered here: the first, referred to as Pr1a, estimates l on the basis of the whole dataset; the second, referred to as Pr1b, estimates l on the basis of the left-censored, considering v 425 m/s so as to meet the condition of independent up-crossing (Table 3). Pr2 is an alternative application that rejects the previous þ assumption of independence. In this case, N V (v) is estimated by

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

31

Fig. 6. Ratio NVþ (v)/fV(v) associated with Grazzanise (a) and probable wind climate (b).

Fig. 7. Gumbel plot of process analysis distributions.

Eq. (8), where l(v) is modeled by an analytical function. The number of simulated wind speeds highlights a regular tendency þ of the ratio NV (v)/fV(v), which has never been detected before due to the lack of sufficiently long recordings in ordinary analyses. To clarify this point further, Fig. 6 shows a plot of the þ þ ratio N V (v)/fV(v), where NV (v) is estimated by counting the upcrossings that have actually occurred, associated with two different datasets: (a) 33 years of 10 min average speeds taken at 3 h intervals at Grazzanise (Italy), i.e. one of the three historical records involved in the definition of the probable wind climate; and (b) 12,738 years of 10 min average speeds at 10 min time intervals, derived from the wind simulation. The ratio in Fig. 6a shows an irregular behavior and its mean value (–) seems to provide a satisfactory estimation. On the contrary, Fig. 6b highlights a regular tendency for v o30 m/s, slight fluctuations around a main trend in the speed range vA[30,37] m/s, and finally unstable behavior for v 437 m/s. The latter is almost definitely due to the approximate estimation of þ NV (v) for high and rare velocities. Assuming that the actual trend of the ratio is captured by the segment associated with vo37 m/ s, and that the terminal unstable behavior is due to the finite size þ of the sample, the ratio N V (v)/fV(v) can be modeled well by the following power law:

Fig. 7 plots the Pr0 (–  –), Pr1 (–) and Pr2 (––) distributions on Gumbel probability paper together with the empirical distribution () of the annual maxima dataset. Comparing Pr0 with the sample distribution, it is possible to detect above what speed range the basic hypothesis of independence between the up-crossings in Eq. (7) is met, that is, v425 m/s. Fig. 7a shows that the bulk of the annual maximum speeds () is accurately fitted by Pr2, while slight detachments can be observed on the tail of the distribution. It is noteworthy that Pr0 and Pr2 are perfectly coincident on the lower tail, thus their departure from the data is ascribable to the unfulfilled hypothesis of independent up-crossings. It seems the data fit very well to the Pr2 distribution in the range of the medium-high velocities. On the contrary, the Pr1a distribution shows a poor fitting. Nevertheless performance improves considerably if the data used for estimating l in Eq. (10) are left-censored, like Pr1b. Fig. 7b highlights that Pr1b (–) is in good agreement with the data in the central body and upper tail of the distribution. On the other hand, the lower tail of Pr1b deviates from the data, even though it remains on the safe side. This behavior of the Pr1b model is ascribable to the previously-observed dependence between up-crossings, which prevents the applicability of process analysis for low speed values.

lðvÞ ¼ dvg

4.3. Storm maxima

ð21Þ

where d and g parameters are reported in Table 3. Fig. 6b shows that l(v) (– –) described by Eq. (21) provides a more accurate þ estimation of NV (v)/fV(v) than its mean value l (–).

The procedure proposed by Cook (1982) to identify independent storms is employed here. A low pass-filter is applied to remove short-duration peaks that may be caused by mechanisms

32

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

parameter FT1 penultimate (– –) distribution also leads to conservative estimates for high return periods, but the amounts of the overestimates are less than those of the FT1 ultimate. Conversely, the use of the single parameter FT1 penultimate distribution (    ) seems to yield design speeds on the unsafe side. The three-parameter FT1 penultimate (–) distribution yields an extremely effective data fitting. Unfortunately, Table 4 shows that, in this case also, the EV distribution resulting from the free fit is not consistent with the parent, as was true during analysis of the annual maxima. Fig. 8a proves that XIMIS (    ) and the three-parameter FT1 penultimate (– –) distributions come to resemble each other since the HW parent is a tail-equivalent Weibull distribution. Hence, if the same regression procedure is applied to estimate the parameters (Table 4), the XIMIS and the three-parameter FT1 penultimate present exactly the same distribution. This means that the XIMIS suffers from the same lack of agreement between parent and extreme distribution as the three-parameter FT1 penultimate model.

other than the dominant one. A lower threshold of 5 m/s is set and each downward crossing of such a threshold is defined as the start of a lull. The resultant average rate of independent storms is 104.7 events per year. This result is very close to the annual rate of independent storms estimated by Cook (1982) for the wind records collected at Kew (104.4) and Lerwick (102.5), while being somewhat lower than that estimated by Cook and Harris (2004) at Boscombe Down (146.5). As observed in a recent paper by Torrielli et al. (2011), the wind climate at Boscombe Down seems to be quite different from the wind climate in the central part of Italy. Similarly to what was done for the annual maxima, only extreme speeds greater than u1 ¼ 22 m/s are fitted to the EV distributions. Such a lower threshold also excludes extreme values associated with secondary mechanisms, which may have got past the low pass-filter. Under these conditions, the analysis of the 386 simulated 33-year speed records leads to the identification of N ¼15,236 wind maxima. Note that analysis of independent storms expands the dataset of annual maxima, which nearly doubles in size. Storm maxima are fitted to the asymptotic FT1 ultimate, FT1 penultimate and XIMIS distributions. The minimum-variance bias-free fitting procedure proposed by Harris (1996) is applied to estimate the parameters of the distributions. Resulting values are listed in Table 4. Fig. 8 presents on a Gumbel plot the fitting of the storm maxima to the analyzed EV distributions. The plotting positions used for the storm maxima () are derived from Eq. (16) setting R¼12,738. FT1 ultimate and penultimate distributions show a better fit to the storm maxima than the annual maxima. However, they show the same tendencies observed in the fitting of the annual maxima. The asymptotic FT1 ultimate (–) tends to be slightly conservative in predicting design speeds with high return periods. The two-

4.4. Other methods using more data from each year Multiple applications of the POT approach are carried out by varying the threshold level u and the fitting procedure used to estimate the GPD parameters, i.e. b and d. The purpose is to draw attention to how those factors affect the probabilistic modeling of the annual maximum speeds. In contrast, the separation time is the same for all applications, i.e. 48 h (Cook, 1985). Fig. 9 shows the distribution of the annual maximum speeds derived using the POT approach. Different values of the threshold u and different regression procedures for estimating the parameters of the GPD distribution are analyzed. In particular, u assumes the values 19, 22, 25, 33 m/s. ML (    ), PWM (– –) and CME (–) regression methods are employed. Table 5 lists the resulting parameters.

Table 4 Parameters of the annual maximum speed distributions from the storm maxima. Values in bold do not result from the regression of the extremes but are assigned in advance. FT1 ultimate

a [m  1  s] u1 [m  s  1] r

0.430 22.44 2559

FT1 penultimate (fit k, a and u1)

FT1 penultimate (fit a and u1)

FT1 penultimate (fit only u1)

XIMIS

k a u1 r

k a u1 r

k a u1 r

w a  ak [m  1  s] P  (au1)k r

1.325 0.196 22.43 1200

1.155 0.277 22.44 3809

1.155 0.324 22.56 20581

Fig. 8. Gumbel plot for annual maximum mean wind speeds from the storm maxima.

1.325 0.115 7.090 1200

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

33

Fig. 9. Distributions of the annual maximum speed derived by using the POT approach. Different thresholds (u) are set.

Table 5 Parameters of the GPD. ML

(a) (b) (c) (d)

PWM

CME

u [m  s  1]

lu [y  1]

b

d [m  s  1]

b

d [m s  1]

b

d [m  s  1]

19 22 25 33

3.761 1.176 0.333 0.009

0.044 0.032 0.023 0.060

2.443 2.320 2.239 2.072

0.072 0.048 0 0

2.508 2.357 2.189 1.956

0.079 0.045 0.036 0.050

2.938 2.587 2.420 2.108

The goodness-of-fit Z-test (Hosking et al., 1984) is used to check the estimate of b derived from the PWM technique. Only in the case (d) the null hypothesis b ¼0 is confirmed. In this case, the high threshold u ¼33 m/s reduces significantly the size of the extremes sample, resulting in a little more than 100 items; this gives rise to strong sampling errors. For medium–low threshold levels, cases (a) and (b), the regression method strongly influences the accuracy of the resulting distribution. On the other hand, for greater thresholds, case (c), different regression methods lead to similar distributions in good agreement with the data. Unfortunately this trend is not confirmed for very high thresholds, case (d), where new discrepancies appear

on both the upper and lower tails of the distributions. These findings seem to suggest that the estimates of the distribution parameters using ML, PWM and CME methods do not converge as the threshold level is raised. Thus, the choice of the fitting procedure and the threshold setting u are highly sensitive issues. Fig. 9 shows that the CME method is the most effective, since it leads to accurate descriptions of the EV distribution in all the 4 cases analyzed. Moreover, good modeling of the annual maxima distributions is obtained by setting u to greater than  25 m/s. In any case, the annual rate of exceedances of the previous thresholds lu is below 1 event/year. This does not meet one of the primary objectives of the POT approach to utilize more data

34

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

for each year than the single annual maximum. The condition lu o1 does not compromise the accuracy of the current analysis, where a large sample of extremes is studied, but it is questionable for real cases. Finally, the present application of the POT method is favored due to the use of speed values in m/s. Real wind data are usually archived in knots, and are affected by rounding errors which may have a great influence on values close to the chosen threshold. As a result, this could penalize the performance of the POT method in real cases. The r-LOS is another EV method that uses more data values for each year. As already stated, r ¼5 can be a suitable number of LOS involved in the EV analysis. It is worth noting that in the current investigation 12,738 years of synthetic wind observations are available (R). It therefore follows that the likelihood function in Eq. (14), or its logarithm, involves more than 63,000 terms, each depending on the three parameters to be estimated. This makes the optimization procedure for evaluating the parameters very burdensome from a computational standpoint. For these reasons, application of the r-LOS method to the present large sized dataset is postponed for further analysis. 4.5. The role of the wind record length The previous analysis, based on a synthetic 12,000-year sample, furnishes a unique item of comparison to evaluate the performance of EV distributions in general terms, but it is not yet able to represent real cases. In this section, the vast numerically-generated database is adopted to assess if the previous conclusions on EV distributions are applicable also when smaller samples of wind observations are available. For this purpose, different cases are explored, characterized respectively by N¼ 25, 50, 100, 200 and 500 years of speed observations. The first two cases represent real situations, where 30–60 years of wind records may be available. The other cases are unrealistic at present, but they are discussed to bridge the real cases with the case N¼ 12,738. For each of the above five cases, a set consisting of 25 records is processed with the aim of evaluating also interval estimates. In this way, not only the expected value (central value of the interval) but also the relative statistical uncertainty (interval width) is given. The EV distributions introduced in the previous paragraphs of Section 4 are used to predict the 1:50 year, 1:200 year and 1:1000 year design speeds: (a) (b) (c) (d) (e)

Asymptotic FT1 ultimate; FT1 penultimate single parameter fitting; XIMIS; GEV; Process analysis by applying both the classical (Pr1b; v425 m/s) and the new (Pr2) forms; (f) POT with u ¼22 m/s.

The FT1 penultimate distribution is applied considering the one-parameter fit, which is the most rigorous form since the parent

is known. The regression procedure stated by Harris is used for estimating the parameters of the FT1 ultimate, penultimate, XIMIS and GEV distributions. In this application of the POT approach, the ML method is used to estimate the parameters of GPD. Both annual and storm maxima are fitted to the FT1 ultimate and penultimate distributions. In light of the findings in the previous paragraphs of Section 4, the set of extreme values are left-censored, and u1 ¼ 25 m/s is assumed to be the lowest allowable speed value. Because of the small size of the dataset involved in the following analysis, significant sampling errors could compromise the choice of the asymptotic form of GEV distribution. To avoid this situation, the Z-test is applied, taking into account also the systematic error (e  0.02). Table 6 reports the frequency with which the GEV distribution results in the FT1, FT2 and FT3 asymptotic forms. Since there are physical grounds for rejecting the FT2 distribution for long-term speed predictions (Holmes and Moriarty, 1999), design speeds deriving from the FT2 distribution are not included in the final interval estimates. As to the XIMIS approach, it is reasonable to expect that the parameter w assumes values within the range [0.5, 3]. Table 6 reports the number of cases in which the estimate of w falls inside and outside such range. In the latter case, the resulting design speeds are rejected from the final interval estimates. The last row of Table 6 concerns the POT approach, and shows the annual number of exceedances over the threshold u ¼22 m/s, averaged over 25 sub-cases of each N. Table 7 reports the 90% confidence interval of the mean of V50, V200 and V1000 derived by applying the previous EV distributions to wind datasets of different sizes (N). It is worth noting that the width of the confidence interval is a measure of the uncertainty with which the central value of the interval is estimated. Thus the width of the interval can be viewed as a measure of the robustness of the EV method. On the one hand, a small width means that the method yields similar predictions from distinct but homogeneous datasets. In this sense the method is considered as robust. The term homogeneous is here referred to wind data associated with the same wind climate. On the other hand, a large width results when different design velocities derive from homogeneous datasets. In this case the method is considered as not robust, since it is sensitive to the differences between distinct but homogeneous datasets. Note that interval estimates are preferred, since the design speed predicted from a single dataset could be not enough representative of the expected values. Fig. 10 plots the 90% confidence interval of the mean of V50 (a) V200 (b) and V1000 (c) against the wind dataset size. The reference speed values (horizontal straight lines) are derived from the empirical distribution relative to the case N¼12,738, i.e. V~ 50 ¼31.23 (a), V~ 200 ¼34.15 (b) and V~ 1000 ¼36.96 m/s (c). Vertical segments define the width of the 90% confidence interval. In general, the 90% confidence intervals become narrower on increasing N. This proves that the statistical uncertainty in prediction of the expected values reduces as the size of the analyzed dataset increases. In any case, the width of the intervals and their location

Table 6 Frequency of occurrence (%) and average values are derived considering a set of 25 items. N¼ 25

N ¼ 50

N ¼100

N ¼ 200

N ¼500

XIMIS

wA[0.5,3] otherwise

28% 72%

52% 48%

60% 40%

72% 28%

92% 8%

GEV

FT1 FT2 FT3

48% 20% 32%

64% 16% 20%

72% 8% 20%

52% 4% 44%

56% 4% 40%

POT

lu [yr  1]

1.457 0.10

1.34 7 0.06

1.25 70.02

1.16 7 0.02

1.17 70.02

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

35

Table 7 90% confidence interval associated with the 1:50, 1:200, 1:1000 yr design values (in m s  1), by analyzing N years of wind values. N¼ 25

N ¼ 50

N ¼100

N ¼ 200

N ¼500

V50 V200 V1000

32.701 7 1.838 36.5687 2.567 41.038 7 3.434

32.4377 1.360 36.1577 1.885 40.458 7 2.500

31.9667 0.797 35.4997 1.086 39.5847 1.427

31.740 70.605 35.185 70.838 39.168 71.109

31.7237 0.401 35.180 7 0.559 39.1767 0.745

V50 V200 V1000

32.1117 1.976 35.6257 2.732 39.6877 3.620

32.019 7 1.146 35.4657 1.565 39.4497 2.055

31.5767 0.726 34.8447 0.990 38.6237 1.300

31.363 70.543 34.539 70.755 38.211 71.002

31.404 7 0.266 34.590 7 0.370 38.2747 0.493

V50 V200 V1000

30.063 7 0.591 32.6687 0.584 35.640 7 0.576

30.0607 0.333 32.6657 0.328 35.6377 0.324

29.9927 0.255 32.5987 0.252 35.5717 0.248

29.961 70.139 32.567 70.138 35.541 70.136

29.9277 0.084 32.5337 0.082 35.507 7 0.081

V50 V200 V1000

30.251 7 0.534 32.8537 0.527 35.8237 0.520

30.276 7 0.295 32.8787 0.291 35.8487 0.287

30.230 7 0.205 32.8327 0.203 35.802 7 0.200

30.220 70.106 32.823 70.105 35.793 70.103

30.237 7 0.071 32.840 7 0.070 35.809 7 0.069

V50 V200 V1000

31.5547 1.349 34.1657 1.768 36.9677 2.280

31.807 7 1.024 34.8677 1.512 38.2927 2.227

31.3377 0.821 34.207 7 1.194 37.3957 1.721

31.250 70.578 34.219 70.900 37.579 71.364

31.206 7 0.337 34.071 7 0.701 37.2877 1.279

V50 V200 V1000

31.2987 1.937 33.7847 3.124 36.4897 4.733

32.054 7 1.524 35.301 7 2.512 38.9317 3.870

31.580 7 0.901 34.660 7 1.549 38.099 7 2.525

31.160 70.780 33.918 71.402 36.884 72.347

31.305 7 0.434 34.2577 0.872 37.4867 1.591

V50

31.4237 1.064

31.4787 0.761

31.330 7 0.572

31.359 70.737

31.2517 0.433

V200 V1000

34.034 7 1.050 37.012 7 1.036

34.088 7 0.751 37.066 7 0.741

33.9427 0.564 36.9217 0.556

33.971 70.728 36.950 70.718

33.8647 0.428 36.8447 0.422

l(v)

V50 V200 V1000

31.1347 0.429 34.020 7 0.440 37.2737 0.452

31.1687 0.287 34.059 7 0.295 37.3177 0.304

31.160 7 0.233 34.053 7 0.242 37.3147 0.252

31.107 70.191 34.001 70.200 37.262 70.208

31.011 7 0.102 33.905 7 0.108 37.1667 0.113

(f) POT

V50 V200 V1000

32.019 7 0.258 34.6657 0.390 37.4417 0.618

32.3987 0.346 35.084 7 0.539 37.8267 0.802

31.6837 0.007 34.4367 0.029 37.3397 0.067

30.985 70.033 33.542 70.068 36.211 70.117

31.3237 0.033 34.1847 0.063 37.2947 0.109

(a) Asymptotic FT1 annual maxima

Storm maxima

(b) FT1 penultimate annual maxima

Storm maxima

(c) XIMIS

(d) GEV (Z-test)

(e) Process analysis

l

relative to the horizontal line, i.e. the reference value, vary from model to model. Below some remarks are reported concerning the EV distributions applied to the synthetic dataset of wind speeds generated in this study: (a) The FT1 ultimate distribution yields over-estimates of the reference values. The amount of the overestimate decreases as N grows, while increasing with the return period from V50 to V1000. The analysis of the storms maxima, instead of the annual maxima, yields expected values (central value of interval) closer to the reference values, and less statistical uncertainties (narrower intervals) for N 425. (b) The FT1 penultimate distribution yields small statistical uncertainties. Unfortunately the resulting design speeds underestimate the reference values. The analysis of annual and storms maxima produces very similar confidence intervals. (c) Generally, the expected value relative to XIMIS predictions is close to the target, or otherwise tends to slight overestimations. This is numerically proved by Table 7. Table 6 explains the reason for the moderate statistical uncertainties. For small data sets (N¼ 25, 50), the shape parameter assumes values outside the reasonable range [0.5, 3] in at least 50% of the cases. The resulting design speed estimates are discarded when evaluating the confidence interval. The statistical uncertainties increase when the size of the analyzed sample diminishes. (d) Similarly to the XIMIS, the estimates from the GEV model also are close to the reference values or sometimes a bit greater. However, the statistical uncertainty is greater than that of the

XIMIS estimates, especially in the prediction of V200 and V1000. This is due to the oscillating form assumed by the GEV model which varies between the FT1 and FT3 asymptotes (Table 6). Cases when the FT2 form is chosen are directly discarded when evaluating the confidence interval so as to reduce the size of the sample analyzed. (e) The Pr1b process analysis yields estimates centered on the reference values. The relative statistical uncertainty is similar to that of the XIMIS distribution in predicting V50, but is less for V200 and V1000. Note that the dispersion of the estimates around the central value is almost independent of the return period. The Pr2 process analysis leads to very close interval estimates, near the reference values. In this case, also, the statistical uncertainty is independent of the return period. (f) The interval estimates from the POT approach are characterized by a very small statistical uncertainty, even though a clear trend cannot be identified. For N ¼100, 200 and 500, the 90% confidence interval collapses for the central value. For small and medium datasets N ¼ 25, 50 and 100, the GPD tends to overestimate the reference values while, on the contrary, slightly underestimating for N ¼200. Predicted values almost coincide with the reference for N ¼500.

5. Conclusions and perspectives This paper proposes a new approach to the hotly debated issue of what is the most efficient EV distribution of annual maximum

36

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

Fig. 10. Prediction of the 1:50, 1:200, 1:1000 yr design values from data set with different size (N).

wind speed for wind engineering applications. A Monte Carlo simulation procedure previously formulated by authors (Torrielli et al. 2010) is adopted to generate artificially unprecedented long-term synthetic wind speed records. They produce 12,738 years of wind speed observations, consistent with the synoptic wind climatology of the central part of Italy. In a first instance, this huge amount of data is fitted to the EV distributions, so as to investigate their performances in modeling annual maximum speed distributions. Next, fixed design speeds are estimated through the empirical distribution of this dataset, and they are used for evaluating the performances of the EV distributions above, when applied to datasets of a size comparable with real situations. These analyses exclude mesoscale wind events for which the EV distributions will be determined separately, and combined with the synoptic EV distributions dealt with in this paper. The application of the unprecedented large sized database provided herein offers the opportunity to draw some remarks on the EV distributions applied to the traditional synthetic dataset of wind observations. Both the analyses of the annual and storms maxima prove that the FT1 ultimate distribution appears to be conservative in predicting design speeds, especially when associated with high return periods. This finding is not supported from a theoretical standpoint. Indeed, in the current study, it is known in advance that the parent distribution is a tail-equivalent Weibull parent,

and so it lies in the domain of attraction of the FT1 asymptote. The highlighted overestimation tendency might be due to the poor convergence of the data to the true asymptotic form, deriving from the forced assumption that rT-N. The FT1 penultimate approach proves to be a robust model, showing a small scatter in the estimates of the design values even in the case of small datasets. Unfortunately, the one-parameter fit FT1 penultimate distribution, which is the most correct in the present application where the parent is completely known, results on the unsafe side. With reference to analysis of large sized datasets, the two-parameter fit FT1 penultimate distribution is quite similar to the FT1 ultimate distribution, sharing the tendency to overestimate design speeds. On the other hand, the three-parameter fit FT1 penultimate distribution leads to a very good data fitting. Nevertheless, the parent parameters estimated from the extreme values are not in agreement with those of the actual parent. Since this issue might be related to the shape parameter of the adopted parent, close to 1, further studies are necessary to investigate this point more deeply. When storms maxima from 12,738-year wind data are used to fit the XIMIS distribution, an accurate description of the annual maxima is obtained. On average, this trend is also for datasets of a size comparable with real cases. However, in 50% of the cases the w shape parameter assumes values outside the expected range [0.5, 3]. Moreover, in the current investigation, the XIMIS approach coincides with the three-parameter fit FT1 penultimate distribution due to the shape of the chosen Hybrid Weibull parent. Therefore, like the parameters of the FT1 penultimate distribution, the parameters of the XIMIS also are not in agreement with the parent parameters. The GEV distribution provides a good description of the annual maxima if it is applied to very large datasets. For smaller datasets, its performances worsen, showing a high statistical uncertainty in predicting design speeds. In addition, the asymptotic form in which the GEV distribution results is not stable over the analyzed datasets. For large datasets, the frequency with which the GEV results in the FT3 asymptote is comparable with that of the FT1 asymptote. These findings contribute to the debate on what is the most correct asymptotic form, with a point scored in favor of the FT3 asymptote. Nevertheless, a series of pending theoretical issues regarding the FT3 asymptote cannot be overlooked, such as the inconsistency with the HW parent, as well as the conditions for the existence of an upper bound of the wind speed (Harris, 2004), an issue not discussed in this paper. The process or level-crossing analysis provides accurate descriptions of the central body and upper tail of the distribution of the annual maximum speeds, but it does not match the description of the lower tail, which has marginal importance for risk analysis. Both classical and new applications of the process analysis yield good results when applied to datasets of 25–50 years of wind data. This approach proves to be robust in predicting design values associated with return periods much longer than the period of observations. It is worth noting that the performance of the classical applications of the process analysis greatly improves if the l parameter is estimated considering only sufficiently high speeds. Besides, the novel use of a power law þ that models the ratio N V (v)/fV(v), in principle very effective, has to be applied extremely carefully after checking a pseudo-regular trend; otherwise it could lead to erroneous results. The POT approach results as being deeply sensitive to decisions relating to the threshold or the method used for parameter evaluation. The CME method has performed well for data fitting, but it is conceptually poor due to the lack of stability of the threshold. Indeed, higher threshold values do not necessarily yield better data fits using the CME method. Besides, a good agreement between regression techniques is obtained for very large datasets and high thresholds. In this case, an effective modeling of the distribution of annual

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

maximum speeds is obtained. Unfortunately, the same results are not obtained in possible situations, when 25–50 years of wind data are available. In these cases the POT approach tends to overestimate like the XIMIS and GEV distributions. In its whole this study points out that some EV methods are more effective than others in the modeling of the distribution of the simulated annual maxima. However, in spite of an accurate data fitting, some theoretical inconsistencies between the EV methods and the parent distribution emerged, to which the authors partially fail to give exhaustive explanations. The authors are aware that some of the findings on the EV distributions obtained in this study are somewhat anomalous with respect to previous analyses but, at the same time, they are convinced of the value of the procedure proposed and of its perspectives. For this reason in this preliminary stage, authors do not express final judgments on the investigated EV distributions, but rather they present to the wind engineering community the outcomes provided by this novel approach, opening such findings to discussion. These outcomes point out the need of further studies currently in progress. They may be grouped into two main classes. The first class of perspectives deals with the need of considering a wider spectrum of wind climatic conditions. It is worth noting that the present analyses are based on synthetic time series of the mean velocity representative of an area in the central part of Italy. This choice was made since the real data from which the analyzed probable wind climate is derived were the best currently available to the authors for performing such a type of analysis. Unfortunately, a wind climate of this type is characterized by a quite unusually low shape parameter of the Weibull parent distribution (k¼ 1.16). In this regard, two lines of development will be pursued. The first concerns a repetition of the analysis with reference to the Ligury region, where k¼1.4–1.7; this aim is particularly straightforward because two of these authors recently created, in the High Tyrrhenian area, a modern monitoring network consisting of 33 ultrasonic anemometers and an acquisition system suitable for creating a huge database of wind speeds recorded continuously with a time step of 0.1 s (Solari et al., 2012); as soon as the quantity of data is large enough, a new wind speed process will be investigated. The second line concerns the attractive opportunity of numerically generating a catalog of wind speed processes coherent with target analytical PDFs and PSDFs, representative of a broad band of worldwide wind climate conditions. The second class of perspectives deals with refining and improving the simulation procedure in order to ensure the generation of mean wind speed time series even closer to real observations. In particular, three goals have to be pursued. The first concerns the parent distribution, in order to establish if the shifting from the Weibull to the Hybrid Weibull model may significantly affect the extreme values drawn from the simulated time series, so as to justify, at least partially, the differences in the EV analyses between this and previous case studies in literature. The second involves the development of a new procedure to take into account the presence of periodic deterministic components in the macro-meteorological wind records (Harris 2008); with such an aim, the simulation procedure has to superimpose periodic deterministic components to the random part of the wind records, while keeping the aim of generating final time series with desired marginal distribution and time correlation. The third concerns the implementation of refined Monte Carlo techniques oriented towards the generation of time series with sample distributions matching almost exactly the target probability distribution; at present such match is almost exact only in terms of ensemble statistics, while slight detaching occurs in the upper tail of the sample distributions. The study of all these topics is in an advanced phase of development and a paper is currently being written to clarify some of these delicate issues.

37

Acknowledgments The present research is part of the European Project ‘‘Wind and Ports’’, financed by the EU Structural Funds of the O.P. Italia– Francia Marittimo 2007–2013. Authors would also acknowledge a reviewer for his great contribution to the improvement of this paper and the development of this research line.

References Abild, J., Andersen, E., Rosbjerg, D., 1992. The climate of extreme winds at the Great Belt, Denmark. Journal of Wind Engineering and Industrial Aerodynamics 41, 521–532. An, Y., Pandey, M.D., 2007. The r-largest order statistics model for extreme wind speed estimation. Journal of Wind Engineering and Industrial Aerodynamics 95 (3), 165–182. Burlando, M., Freda, A., Ratto, C.F., Solari, G., 2010. A pilot study of the wind speed along the Rome–Naples HS/HC railway line. Part 1-Numerical modelling and wind simulations. Journal of Wind Engineering and Industrial Aerodynamics 98 (8-9), 392–403. Cauchy, A.L.,1821. Cours d’analyse I. Coles, S.G., 2001. An introduction to Statistical Modelling of Extreme Values. Springer, Berlin. Coles, S.G., Walshaw, D., 1994. Directional modelling of extreme wind speeds. Applied Statistics 43, 139–157. Cook, N.J., 1982. Towards better estimation of extreme winds. Journal of Wind Engineering and Industrial Aerodynamics 9, 295–323. Cook, N.J., 1985. The Designer’s Guide to Wind Loading of Building Structures. Part 1: Background, Damage Survey, Wind Data and Structural Classification. Building. Research Establishment, Garston and Butterworths, London 371. Cook, N.J., 2004. Confidence limits for extreme wind speeds in mixed climates. Journal of Wind Engineering and Industrial Aerodynamics 92, 41–51. Cook, J., Harris, I., 2004. Exact and general FT1 penultimate distributions of extreme wind speeds drawn from tail-equivalent Weibull parents. Structural Safety 26, 391–420. Cook, J., Harris, I., 2008. Postscript to ‘‘Exact and general FT1 penultimate distributions of extreme wind speeds drawn from tail-equivalent Weibull parents’’. Structural Safety 30, 1–10. Cramer, H., 1946. Mathematical Methods of Statistics. Princeton University Press, Princeton NJ. Davenport, A.G., 1968. The dependence of wind loads on meteorological parameters. In: Proceedings of the International Conference on Wind Effects on Buildings and Structures, Ottawa 1, pp. 19–82. Davison, A.C., 1984. Modelling excess over high thresholds with an application. In: de Oliveira, J.Tiago (Ed.), Statistical Extremes and Applications. Reidel Dordrecht, pp. 461–482. Davison, A.C., Smith, R.L., 1990. Models for exceedances over high thresholds. Journal of the Royal Statistical Society B 52, 393–442. Fisher, R.A., Tippett, L.H.C., 1928. Limiting forms of the frequency distribution of the largest or smallest member of a sample. Proceedings of the Cambridge Philosophical Society 24, 180–190. Freda, A., Solari, G., 2010. A pilot study of the wind speed along the Rome–Naples HS/HC railway line. Part 2 – Probabilistic analyses and methodology assessment. Journal of Wind Engineering and Industrial Aerodynamics 98, 404–416. Gomes, L., Vickery, B.J., 1977. On the prediction of extreme wind speeds from the parent distribution. Journal of Wind Engineering and Industrial Aerodynamic 2, 21–36. Grigoriu, M., 1995. Applied Non-Gaussian Processes: Examples, Theory, Simulation, Linear Random Vibration, and Matlab Solutions. Prentice-Hall, Englewood Cliffs, NJ. Gumbel, E.J., 1958. Statistics of Extreme. Columbia University, New York. Harris, R.I., 1996. Gumbel re-visited - a new look at extreme value statistics applied to wind speeds. Journal of Wind Engineering and Industrial Aerodynamic 59 (1), 1–22. Harris, R.I., 1999. Improvements to the ‘Method of Independent Storms’. Journal of Wind Engineering and Industrial Aerodynamic 80, 1–30. Harris, R.I., 2004. Extreme value analysis of epoch maxima—convergence, and choice of asymptote. Journal of Wind Engineering and Industrial Aerodynamic 92, 341–360. Harris, R.I., 2005. Generalized Pareto methods for wind extremes. Useful tool or mathematical mirage? Journal of Wind Engineering and Industrial Aerodynamic 93, 897–918. Harris, R.I., 2006. Errors in GEV analysis of wind epoch maxima from Weibull parents. Wind and Structures 9 (3), 171–191. Harris, R.I., 2008. The macro-meteorological spectrum—a preliminary study. Journal of Wind Engineering And Industrial Aerodynamics 96, 2294–2307. Harris, R.I., 2009. XIMIS, a penultimate extreme value method suitable for all types of wind climate. Journal of Wind Engineering and Industrial Aerodynamic 97, 271–286. Holmes, J.D., 2003. Discussion of ‘Generalized extreme gust wind speeds distributions’ by E. Cheng and C. Yeung. Journal of Wind Engineering and Industrial Aerodynamic 91, 956–957.

38

A. Torrielli et al. / J. Wind Eng. Ind. Aerodyn. 115 (2013) 22–38

Holmes, J.D., Moriarty, W.W., 1999. Application of the generalized Pareto distribution to extreme value analysis in wind engineering. Journal of Wind Engineering and Industrial Aerodynamic 83, 1–10. Hosking, J.R.M., Wallis, J.R., Wood, E.F., 1984. Estimation of the generalized extreme-value distribution by the method of probability weighted moment. Institute of Hydrology, Report, 89. Jenkinson, A.F., 1955. The frequency distribution of the annual maximum (or minimum) values of meteorological elements. Quarterly Journal of the Royal Meteorological Society, 86. Kai, K., 1987. Spectrum climatology of the surface winds in Japan. Part II: the diurnal variations, the synoptic fluctuations, the 10–20 days fluctuations and the annual variation. Journal of the Meteorology Society of Japan 65 (5), 751–765. Lieblein, J., 1974. Efficient method of extreme-values methodology. National Bureau of Standards, Washington DC, Report NBSIR 74-602, 31. Masters, F., Gurley, K.R., 2003. Non-Gaussian simulation: cumulative distribution function map-based spectral correction. Journal of Engineering Mechanics 129 (12), 1418–1428. Naess, A., 1998. Estimation of Long Return Period Design Values for Wind Speeds. Journal of Wind Engineering and Industrial Aerodynamic 124 (3), 252–259. Palutikof, J.P., Brabson, B.B., Lister, D.H., Adcock, S.T., 1999. A review of methods to calculate extreme wind speeds. Meteorological Applications 6, 119–132. Peterka, J.A., 1992. Improved extreme wind prediction for the United States. Journal of Wind Engineering and Industrial Aerodynamic 41-44, 533–541. Pickands, J., 1975. Statistical inference using extreme order statistics. Annals of Statistics 3, 119–131. Rice, S.O., 1944. Mathematical analysis of a random noise. Bell System Technical Journal 23, 282–332. Rice, S.O., 1945. Mathematical analysis of a random noise. Bell System Technical Journal 24, 46–156. Simiu, E., Heckert, N.A., 1996. Extreme wind distribution tails - a peaks over threshold approach. Journal of Structural Engineering ASCE 122, 539–547. Simiu, E., Filliben, J.J., 1999. Discussion of ‘‘Design Gust Wind Speeds in the United States. Journal of Structural Engineering ASCE 125 (9), 1077–1078. Smith, R.L., 1986. Extreme value theory based on the r largest annual events. Journal of Hydrology 86, 27–43.

Solari, G., Repetto, M.P., Burlando, M., De Gaetano, P., Pizzo, M., Tizzi, M., Parodi, M., 2012. The wind forecast for safety and management of port areas. Journal of Wind Engineering and Industrial Aerodynamic 104-106, 266–277. Takle, E.S., Brown, J.M., 1978. Note on the use of Weibull statistics to characterize wind speed data. Journal Of Applied Meteorology 17, 556–559. Tawn, J.A., 1988. An extreme-value theory model for dependent observations. Journal of Hydrology 101, 227–250. Torrielli, A., Repetto, M.P., Solari, G., 2010. Simulation and analysis of long-period samples of mean wind velocity. In: Proceedings of the 5th International Symposium on Computational Wind Engineering, Chapel Hill, North Carolina. Torrielli, A., Repetto, M.P., Solari, G., 2011. Long-term simulations of the mean wind velocity. Journal of Wind Engineering and Industrial Aerodynamic 99, 1139–1150. Torrielli, A., 2011. Long term simulation and reliability analysis of structures subjected to the wind action. Ph.D. Thesis. Available from: /http://www. scuoladottorato-sti.unige.it/content/view/134/156/lang,italian/S. Torrielli, A., Repetto, M.P., Solari, G. A refined simulation for investigating extreme wind speeds. Journal of Wind Engineering and Industrial Aerodynamic, in preparation. Von Mises, R., 1936. La distribution de la plus grandes de n valeurs. Reprinted in Selected Papers II, American Mathematical Society, Providence, RI (1954), pp. 271–294. Walshaw, D., 1994. Getting the most from your extreme wind data: a step by step guide. Journal of research of the National Institute of Standards and Technology 99, 399–411. Weiss, L., 1971. Asymptotic inference about a density function at the end of its range. Naval Research Logistic Quarterly 18, 111–114. Weissman, I., 1978. Estimation of parameters and large quantities based on the k largest observations. Journal of the American Statistical Association 73, 812–815. Young, D.H., 1967. Recurrence relations between the PDF’s of order statistics of dependent variables, and some applications. Biometrika 54 (1/2), 283–292. Yu, G.-H., Huang, C.-C., 2001. A distribution free plotting position. Stochastic Environmental Research and Risk Assessment 15, 462–476.