Simulation of non-Gaussian field applied to wind pressure fluctuations

Probabilistic Engineering Mechanics 15 (2000) 339–345 www.elsevier.com/locate/probengmech

Simulation of non-Gaussian field applied to wind pressure fluctuations M. Gioffre` a,*, V. Gusella a, M. Grigoriu b a

Department of Water and Structural Engineering, University of Perugia, Via G.Duranti 93, 06125 Perugia, Italy b School of Civil and Environmental Engineering, Cornell University, Ithaca, NY 14853, USA

Abstract A simulation algorithm to generate non-Gaussian wind pressure fields is proposed. This algorithm uses the correlation–distortion method based on translation vector processes. Conditions on the matrix of cross-covariance functions are given to assure the applicability of the model. The proposed method does not require iterative procedures and it is well suited when experimental data are available. In particular it requires cross-covariance functions and marginal distribution that can be directly estimated from data. To illustrate the procedure, the model is calibrated on experimental results obtained from wind tunnel tests on a tall building. The efficiency of the proposed methodology for reproducing the non-Gaussian nature of pressure fluctuations on separated flow regions is demonstrated. q 2000 Elsevier Science Ltd. All rights reserved. Keywords: Wind pressure field; Non-Gaussian vector process; Translation vector process; Non-Gaussian field simulation

1. Introduction Wind induced aerodynamic loads are strictly connected to the flow-structure interactions. Several studies have already outlined that the Gaussian model can be used with good accuracy in modeling the pressure fluctuations in windward regions of buildings where the flow is attached. On the other hand, the pressure field clearly shows its nonGaussian features in the regions that follows a separation of the flow [1–3]. This separation is influenced by several parameters like the incoming wind direction and the superficial roughness of the body. It is well known that the nonGaussian pressure field in the separated flow regions of the bluff bodies, which are characteristic structures of interest in civil engineering, cannot be adequately modeled yet by theoretical and/or numerical models. Furthermore, the accurate prediction of all the aerodynamic loads is very important for estimating the loads and load effects produced by extreme winds as well as the fatigue damage of claddings and fastening. To overcome the difficulties associated with the aerodynamic loads prediction, the combination of experimental tests and simulation procedures can be a very useful tool. This approach has been followed by several researchers in the last years [4,5]. Nevertheless, the accuracy and efficiency in the generation of non-Gaussian * Corresponding author. Tel.: 139-075-585-2616; fax: 139-075-5852830. E-mail addresses: [email protected] (M. Gioffre`), [email protected] (V. Gusella).

wind pressure time histories is still being pursued. Within this context, an algorithm for generating multivariate nonGaussian processes is proposed in this paper. The procedure uses the correlation–distortion method based on translation vector processes. The proposed algorithm is applied to the generation of wind pressure time series in the separated flow regions of a prismatic tall building. The experimental data used to calibrate the model are constituted of the wind pressure fluctuations recorded in wind tunnel tests [6]. 2. Translation vector process The wind pressure fluctuations acting on the surfaces of a building can be modeled by a stochastic time dependent field P ˆ P…z1 ; z2 ; z3 ; t†

…1†

where z ˆ {z1 ; z2 ; z3 }T is the vector with components, which represent the coordinates of the considered point in space, and t is the time. From Eq. (1) it is possible to introduce the coefficient of pressure field Cp ˆ Cp …z1 ; z2 ; z3 ; t† removing the static component of the pressure fluctuations and normalizing with the dynamic pressure at a reference height. The coefficient of pressure fluctuations Cp …z; t† can be modeled by a one-variate four-dimensional (1V, 4D) stochastic non-Gaussian field. Information on the pressure field is generally available at a discrete number of points, e.g. the m tap locations where the pressure is measured in experimental tests. The collection of the coefficient of

0266-8920/00/$ - see front matter q 2000 Elsevier Science Ltd. All rights reserved. PII: S0266-892 0(99)00035-1

M. Gioffre` et al. / Probabilistic Engineering Mechanics 15 (2000) 339–345

340

pressure fluctuations Cp1 …t†; …; Cpm …t† into a vector allows the description of the pressure fluctuations by a discrete stochastic field represented by an m-variate one dimensional (mV, 1D) stochastic vector process. Thus, with this discretization, the 1V, 4D stochastic field is described by a mV, 1D stochastic vector process. The non-Gaussian nature of the wind pressure vector process on tall buildings can be successfully described by a translation model. Let Y…t† ˆ {Y1 …t†; …; Ym …t†}T ; t [ R; be an ergodic m-variate Gaussian vector process with zero vector mean, covariance functions rkl …t† ˆ E‰Yk …t†Yl …t 1 t†Š; k; l ˆ 1; …; m; and marginal distributions equal to the standard Gaussian distribution F Yk …yk † ˆ P‰Yk …t† # yk Š ˆ F Y …yk †; k ˆ 1; …; m: The translation vector process X…t† ˆ {X1 …t†; …; Xm …t†} T is defined by the non-linear memoryless transformation g : Rm ! Rm with monotonic components. The kth component of X(t) is given by [7] ‰F Y …Yk …t††Š Xk …t† ˆ gk ‰Yk …t†Š ˆ FX21 k

…2†

where FXk …xk † is the marginal distribution of the kth component Xk(t) of X(t). The vector translation process X(t) is completely characterized by the mapping g(·) and the covariance functions rkl …t†; k; l ˆ 1; …; m; because Y(t) is a Gaussian vector process. The vector mean mX, and covariance functions CXk Xl …t† ˆ E‰Xk …t†Xl …t 1 t†Š 2 mXk mXl of X(t) can be directly calculated using the definitions in Eq. (2). The scaled covariance functions of X(t) are defined as

j kl …t† ˆ

CXk Xl …t† sX k sX l

…3†

k; l ˆ 1; …; m; and represent the covariance functions of the vector process with components X~ k …t† ˆ ‰Xk …t† 2 mXk Š=sXk : 2.1. Properties of the covariance functions Arbitrarily defined scaled covariance functions and marginal distributions are not always compatible with a description of the non-Gaussian vector process X(t) by means of a translation model. Let fYk1 Yl2 ‰yk1 ; yl2 ; rkl …t†Š be the joint density of the dependent Gaussian variables Yk(t) and Yl …t 1 t†: The covariance functions j kl(t ) and r kl(t ) are related by the equation Z∞ Z∞ s Xk sXl j kl …t† ˆ ‰gk …yk1 † 2 mxk Š‰gl …yl2 † 2∞

2∞

2 mxl ŠfYk1 Yl2 ‰yk1 ; yl2 ; rkl …t†Š dyk1 dyl2

…4†

k; l ˆ 1; …; m Eq. (4) can be used to calculate the values of j kl(t ) corresponding to r kl(t ) for an arbitrary but fixed t . Furthermore, the following inequality holds [8]: uj kl …t†u # urkl …t†u

k; l ˆ 1; …; m

…5†

Considering the auto-covariance functions one has that, if

rkk …t†; k ˆ 1; …; m; is equal to 1, 0, or 21, the corresponding values of j kk(t ) are 1, 0, or the “lower bound” [7] pl j kk ˆ

E‰gk …Yk †gk …2Yk †Š 2 {E‰gk …Yk †Š}2 E‰gk …Yk †2 Š 2 {E‰gk …Yk †Š}2

…6†

respectively, where Yk …t† ˆ 2Yk …t 1 t† ˆ Yk : If one considers the cross-covariance functions (k ± l) it happens that even the limit rkl …t† ˆ 1 does not always imply that j kl …t† ˆ 1: For example, assume the two-dimensional Gaussian vector with components Y1 ˆ Y2 ˆ Y with mean zero and unit variance. Define the translation vector components as X1 ˆ g1 …Y1 † ˆ Y and X2 ˆ g2 …Y2 † ˆ Y 2 : The variables Y1 and Y2 are perfectly correlated but the expectation E‰X1 …X2 2 1†Š is equal to E‰Y 3 Š 2 E‰YŠ ˆ 0: Let FX0 k …xk † and j kl0 …t†; k; l ˆ 1; …; m; be the target marginal distributions and covariance functions. The covariance functions r0kl …t† have to be found using Eq. (4). This equation has solution if every covariance j kl0 …t† is in the interval ‰j klpl ; j klpu Š for every t . The lower bound j klpl and the upper bound j klpu are calculated corresponding to rkl …t† ˆ 21 and rkl …t† ˆ 1; respectively. If solutions r0kl …t† to Eq. (4) exist, one has to verify the functions r0kk …t† to be non-negative definite, i.e. aT r…n† 0kk a $ 0

k ˆ 1; …; m

…7†

where n is a positive integer, a [ Rn denotes an arbitrary 0 vector, and r…n† 0kk ˆ {rkk …ti 2 tj †}; i; j ˆ 1; …; n; denotes a matrix consisting of the covariances of the components Yk …ti † and Yk …tj † of the random vector Yk ˆ {Yk …t1 †; …; Yk …tn †}: Alternative ways of checking whether or not r…n† 0kk is non-negative definite are

lmin ‰r…n† 0kk Š $ 0

…8a†

ur…n† 0kk u $ 0

…8b†

where

lmin ‰r…n† 0kk Š

and ur…n† 0kk u denote the smallest determinant of r…n† 0kk ; respectively.

eigenvalue

and the Furthermore, it has to be verified if the matrix 1 0 0 r11 …t† … r01m …t† C B B . .. C C . r0 … t † ˆ B B . ] . C A @

…9†

r0m1 …t† … r0mm …t† is non-negative definite for every time lag t . If at least one of the r0kk …t† and/or r0 …t† are negative definite, then they are not covariances and it is not possible to represent the target non-Gaussian process as a translation vector process with scaled covariances j kl0 …t† and marginal distributions FX0 k …xk †: 3. Multivariate simulation The objective of the proposed procedure is the generation of samples of the m-variate non-Gaussian vector process X(t) with the target matrix of the m × m scaled covariance

M. Gioffre` et al. / Probabilistic Engineering Mechanics 15 (2000) 339–345

representation is given by

functions 0 B B j …t† ˆ B B @ 0

341

0 j 11 …t†

.. .

… ]

0 j 1m …t†

.. .

1

Yk;wk …t† ˆ

C C C C A

…10†

wk;t 1 wk 1 1 X uˆ2 wk;t 2 wk

Yk;u au …t; Tk †

…11†

wk;t Tk # t # …wk;t 1 1†Tk

0 0 j m1 …t† … j mm …t†

FX0 k …xk †;

and the target marginal distributions k ˆ 1; …; m of the kth component of X(t). The procedure has three main steps: model calibration, sample generation, and mapping. The model calibration consists in determining: • the second-moment characteristics of the Gaussian image Y(t) of X(t); • the m memoryless non-linear transformations gk …·†; such that the target marginal covariance matrix j 0(t ) and marginal distributions F X0k …xk † can be matched. Taking FXk …xk † ˆ FX0 k …xk † in Eq. (2) the vector non-linear transformation g(·) is completely determined and the target marginal distributions are matched. It is worthwhile to note that this method is particularly attractive for pressure simulation since it is possible to match exactly the marginal distributions which can be very different from one tap location to another depending on the complexity of flow– structure interactions. Once the transformation g(·) is set, the conditions on the covariances have to be verified. First, the upper and lower bounds, j klpu and j klpl ; k; l ˆ 1; …; m; have to be determined. If j kl0 …t† [ ‰j klpl ; j klpu Š; ;t then the functions r0kl …t† can be calculated, analytically or numerically, using Eq. (4). The successive verification is on the non-negative definiteness of the functions r0kk …t†; k ˆ 1; …; m; and the matrix r 0(t ). The non-Gaussian vector process X(t) can be modeled by a vector translation process, in matching j 0(t ) and F X0k …xk †; only if all these conditions are verified. Otherwise, two main strategies can be adopted. The first one is considering more general classes of non-Gaussian models that can better represent the desired statistics [7]. The second is approximating j 0(t ) and/or F X0k …xk †: The particular approximation to be considered depends on the nature of the problem. For example, if one is more interested in calculating extremes, j 0(t ) can be approximated and F X0k …xk † can be matched exactly. Once all the conditions are verified, the second step of the simulation algorithm consist of generating realizations of the m-variate Gaussian process Y(t) with cross-covariance matrix r 0(t ). Either one of the well established methods reported in literature can be used. In the present work a method based on the sampling theorem has been preferred [1] because it appears to be the most efficient since it only requires the knowledge of the covariance functions at discrete times, which can be directly estimated from data. This method utilizes a local parametric representation of the components Yk …t† of Y(t), with frequency bands ‰0; vpk Š: The

k ˆ 1; …; m; where Tk ˆ p=vpk ; wk is a positive integer determining the size of the window used in the representation, wk;t is the largest integer smaller than t=Tk ; Yk;u ˆ Yk …uTk † are values of component k of Y(t) at equally spaced times uTk, called nodal points, and

au …t; Tk † ˆ

sin‰p…t 2 uTk †=Tk Š p…t 2 uTk †=Tk

…12†

Background on this model is beyond the scope of this paper and can be found in [9]. Samples of Y(t) are generated using the numerical procedure based on this model and proposed in the work of Gioffre` [1]. In the third and last step of the algorithm the generated samples of Y(t) are transformed into samples of X(t). Each component yk …t† is mapped to the corresponding xk …t† using the set of equations given by Eq. (2). The main characteristics of the proposed procedure can be summarized as follows: 1. Exact determination of r 0(t ) by numerical or analytical tools once its existence is verified; 2. No need of iterative procedures because r 0(t ) is calculated only once prior to the generation of the Gaussian samples; 3. The efficiency of the algorithm is only related on the effectiveness of the Gaussian simulation. It is worth noting that the first two characteristics represent also the main advantages as opposed to other methods based on the correlation–distortion method reported in literature. In particular the methods reported in [10,11] use approximate relations between j kl …t† and rkl …t† (or, equivalently, between their Fourier transforms) based on Hermite polynomials and a limited number of statistical moments. Furthermore, the methods reported in [12,13] require several iterative generation of Y(t) to match the target correlation functions. 4. Wind pressure vector field 4.1. Experimental data Fig. 1 shows a view of the CAARC standard tall building model utilized in the wind tunnel experimental tests performed in the boundary layer wind tunnel of CRIACIV (Italian Inter-university Research Center in Aerodynamics of Constructions and Wind Engineering) located in Prato, Italy [1]. The data used for this work have been obtained using the simple model experiment proposed by the Commonwealth Advisory Aeronautical Research Council

342

M. Gioffre` et al. / Probabilistic Engineering Mechanics 15 (2000) 339–345

Fig. 2. Cross-section of the model. Taps numbering at 2/3h.

The coefficient of pressure is defined by Cp …t† ˆ

p…t† 2 ph …t† …1=2†ra U…h† 2

…13†

were p(t) is the total pressure, h the building model height, ph …t† the static pressure at h, r a the mass density of air, and U(h) the mean wind velocity at h. The three Cp(t) time series were recorded in the separated flow region at taps 26, 28 and 30 located at 2/3 the height h of the model, with incoming wind normal to the north face (b ˆ 08 in Fig. 2). Plots of the recorded time histories are shown in Figs. 3–5 where the non-Gaussian features are evident.

Fig. 1. View of the CAARC model and pressure taps numbering.

(CAARC) in 1969. According to CAARC specifications [14] the boundary layer was to be representative of wind blowing over a long fetch of forest or urban development with building heights in the range of 6–15 m. The geometry of the model was prismatic with a constant rectangular cross-sections and fixed ratios between dimensions of 1:1.5:6. The scale used for the model was 1:500 and the dimensions were dx ˆ 61 mm; dy ˆ 91 mm; and h ˆ 366 mm: Sixty tap locations were placed on the model and the wind pressures were recorded for different wind directions (angle b in Fig. 1). The frequency response of the pressure measurement system was essentially flat up to 100 Hz. Three time series of coefficient of pressure fluctuations, Cp …t†; are used in this paper to apply the proposed procedure.

Fig. 3. Recorded experimental time series at tap 26, simulated Gaussian sample y1(t) and corresponding non-Gaussian sample xp1 …t† after mapping.

M. Gioffre` et al. / Probabilistic Engineering Mechanics 15 (2000) 339–345

343

Table 1 Mean, standard deviation, skewness and kurtosis coefficients

X1p …t† X1(t) XL1(t) X2p …t† X2(t) XL2(t) X3p …t† X3(t) XL3(t)

m

s

g3

g4

20.6637 0 0 20.7285 0 0 20.7371 0 0

0.2828 1 1 0.3313 1 1 0.3937 1 1

21.1281 1.1281 1.1281 21.2736 1.2736 1.2736 21.9572 1.9572 1.9572

5.0811 5.0811 5.3452 6.0423 6.0423 6.0162 4.5396 4.5396 4.6724

the marginal distributions F X0k …xk †; k ˆ 1; 2; 3 of the translation model. The covariance functions are shown in Fig. 8 (continuous line). The shifted Lognormal translation processes Fig. 4. Recorded experimental time series at tap 28, simulated Gaussian sample y2(t) and corresponding non-Gaussian sample xp2 …t† after mapping.

XLk …t† ˆ ak 1 exp‰mYkp 1 sYkp Y…t†Š

k ˆ 1; 2; 3

…14†

The Cp;j …t†; j ˆ 26; 28; 30 are assumed to be the components Xkp …t†; k ˆ 1; 2; 3; of an ergodic vector process X p(t). The mean mXkp ; standard deviation sXkp ; skewness coefficient gp3;k ; and kurtosis coefficient gp4;k ; of the components Xkp …t† are reported in Table 1. Let X~ k …t† ˆ …Xkp …t† 2 mXkp †=sXkp ; k ˆ 1; 2; 3; be the three ~ standardized components of X…t†: For simplicity of notation ~ the three-dimensional vector process X…t† ˆ 2X…t† is considered in the following. The normalized time series Xk …t† have been used to estimate the matrix of target scaled covariances j 0(t ) and

where ak ; mYkp ; and sYkp are the parameters of the model, have been found to be the better fit with the assumed time series. The first four moments of the Lognormal distributions, assumed as the non-Gaussian model, are reported in Table 1. The shape of the three marginal densities is shown in Figs. 6 and 7 (continuous lines). Given the non-linear transformations gk …Y…t†† of Eq. (14), it has been possible to calculate the lower and upper bounds of the covariances. It was found that j kl0 …t† [ ‰j klpl ; j klpu Š; k; l ˆ 1; …; m; for every t . It follows that a unique solution r 0(t ) exists and can be calculated using Eq. (4). The functions r0kk …t†; k ˆ 1; 2; 3; and the complete matrix r 0(t ) have been found to be positive definite. All the required conditions were therefore satisfied and it was

Fig. 5. Recorded experimental time series at tap 30, simulated Gaussian sample y3(t) and corresponding non-Gaussian sample xp3 …t† after mapping.

Fig. 6. Estimated marginal densities from simulation and marginal densities of the Lognormal model.

4.2. Model calibration

344

M. Gioffre` et al. / Probabilistic Engineering Mechanics 15 (2000) 339–345

possible to represent the vector process X(t) by means of a translation model to match exactly the target marginal distributions FX0 k …xk †; k ˆ 1; 2; 3; and covariance matrix j 0 …t†:

5. Numerical results

Fig. 7. Estimated marginal densities from simulation and marginal densities of the Lognormal model in logarithmic scale.

Samples y(t) of the 3-variate Gaussian process Y(t) with zero mean and covariance matrix r 0(t ) are generated using the method based on the parametric representation in Eq. (11). It is assumed that w k ˆ w ˆ 3; and Tk ˆ T ˆ 0:004 s, k ˆ 1; 2; 3: The obtained samples are mapped into samples x(t) of X(t) by the transformations (14). Figs. 3–5 show realizations of the components of the Gaussian vector Y(t) and the Lognormal vector X p(t) after mapping. The comparison between the behavior of the simulated time series and the recorded data from which the model was calibrated is satisfactory. Fig. 6 compares the marginal density functions of the XLk ; k ˆ 1; 2; 3 (continuous lines) and the estimation of the marginal histograms (W) of the components Xk …t† of X(t) based on 100,000 realizations. Plotting the densities

Fig. 8. Target (continuous lines) and estimated from simulation (W) covariance functions of the standardized time series X~ k …t†; k ˆ 1; 2; 3:

M. Gioffre` et al. / Probabilistic Engineering Mechanics 15 (2000) 339–345

of each single tap in logarithmic scale it is evident that the estimated density from simulation is in agreement with the model also in the tail regions (Fig. 7). In Fig. 8 the circles represent the estimated covariance functions based on 10,000 realizations of the vector X(t), which perfectly match j 0(t ) (continuous lines). 6. Conclusions A procedure for generating multivariate wind pressure time histories has been presented. The procedure uses the correlation–distortion method based on translation vector processes. The mathematical background and the conditions for the applicability of the algorithm have been reported. The proposed procedure has been applied to the generation of the wind pressure time series in the separated flow region of a tall building. It has been verified that a Lognormal translation vector process accurately describes the main non-Guassian features of the wind pressure field. Moreover, the conditions for the applicability of the proposed procedure are satisfied in this case study. It was therefore possible to generate samples of wind pressure time series avoiding time consuming iterative procedures. Furthermore, if the required conditions are not verified, it is always possible to control the accuracy in the description of the second-moment characteristics and/or the marginal distribution of the target nonGaussian vector. Depending on the particular problem study considered one can decide to approximate one or both j 0(t ) and F X0k …xk †; k ˆ 1; …; m: It has been shown that the efficiency of the proposed procedure is strictly related to the simulation of the parent Gaussian vector process. In this paper a method based on the sampling theorem is used because it only requires the knowledge of the covariance function at discrete times, which is directly estimated from experimental data. The analysis of the marginal distributions and the crosscovariance functions, estimated from the simulation of three time histories at the vortex shedding region of a tall building,

345

demonstrated that the accuracy and the efficiency of the proposed procedure are very satisfactory. References [1] Gioffre` M. Analysis and simulation of non-Gaussian processes with application to wind engineering and reliability. PhD thesis, University of Florence, Florence, Italy, 1998. [2] Peterka JA, Cermak JE. Wind pressure on buildings—probability densities. Journal of the Structural Division 1975;101(ST6):1255–67. [3] Stathopoulos T. Pdf of wind pressures on low-rise buildings. Journal of the Structural Division 1980;106(ST5):973–90. [4] Gurley KR, Tognarelli MA, Kareem A. Analysis and simulation tools for wind engineering. Probabilistic Engineering Mechanics 1997;12(1):9–31. [5] Kumar KS, Stathopoulos T. Computer simulation of fluctuating wind pressures on low building roofs. Journal of Wind Engineering and Industrial Aerodynamics 1997;71:485–95. [6] Gioffre` M, Gusella V, Grigoriu M. Analysis of non-Gaussian stochastic field of wind pressure. In: Proceedings, Fourth International Conference on Stochastic Structural Dynamics, Notre Dame, IN, 1998. [7] Grigoriu M. Applied non-Gaussian processes: examples, theory, simulation, linear random vibration and MATLAB solutions. Englewood Cliffs, NJ: Prentice-Hall, 1995. [8] Lancaster HO. Some properties of the bivariate normal distribution considered in the form of a contingency table. Biometrika 1957;44:289–92. [9] Grigoriu M, Balopoulou S. A simulation method for stationary Gaussian random functions based on the sampling theorem. Probabilistic Engineering Mechanics 1993;8(3–4):239–54. [10] Gurley KR, Kareem A. Simulation of correlated non-Gaussian pressure fields. MECCANICA 1998;33(3):309–17. [11] Gurley KR, Kareem A, Tognarelli MA. Simulation of a class of nonnormal random processes. Journal of Non-Linear Mechanics 1996;31(5):601–17. [12] Popescu R, Deodatis G, Prevost JH. Simulation of homogeneous nonGaussian stochasti vector fields. Probabilistic Engineering Mechanics 1998;13(1):1–13. [13] Yamazaki F, Shinozuka M. Digital generation of non-Gaussian stochastic fields. Journal of Engineering Mechanics, ASCE 1988;114(7):1183–97. [14] Wardlaw RL, Moss GF. A standard tall building model for the comparison of simulated natural winds in wind tunnels. Report CC662 Tech. 25, Commonwealth Advisory Aeronautical Research Council, January 1970.