Numerical study of pressure fluctuations on a rectangular cylinder in aerodynamic oscillation

ELSEVIER

Journal of Wind Engineering and Industrial Aerodynamics54/55 (1995) 239-250

~ ~ g ~

Numerical study of pressure fluctuations on a rectangular cylinder in aerodynamic oscillation T. T a m u r a a, Y. I t o h a, A. W a d a ~, K. K u w a h a r a b a Tokyo Institute of Technology, 4259 Nagatsuta, Midori-ku, Yokohama 227, Japan b The Institute of Space and Astronautical Science, Kanagawa 229, Japan

Abstract

We numerically simulate unsteady flow fields around a forced-oscillating rectangular cylinder. The higher-order finite difference computation of the two-dimensional incompressible Navier-Stokes equations is applied to the dynamic interaction between the vortex motions and the torsional oscillation of a rectangular cylinder. In order to clarify the physical mechanism of the aerodynamic damping forces and the self-excited forces, we investigate the pressure fluctuations acting on the side surface. According to the computational flow visualization of the instantaneous vorticity, the vortices shed from the leading edge are controlled by the motion of a rectangular cylinder and the lock-on patterns of vortex shedding are dearly seen under the proper conditions.

1. Introduction

The unstable oscillation of rectangular cylinders occurs owing to the self-excited force mainly stemming from interactions with the separated shear flows or vortices close to the cylinder. The flow pattern with separation region around an oscillating rectangular cylinder is so complicated that the analytical theory of aeroelasticity cannot be easily applied to this bluff shape of the section and the experimental approach to this problem is also difficult for comprehending the physical mechanism of the aerodynamic instability due to the separated shear flows. Thus so far we have computed three-dimensional incompressible flows around a square cylinder at rest and its aerodynamic characteristics and the importance of three-dimensional simulations has been clarified [1,2]. We also simulated two-dimensional unsteady flows around an oscillating rectangular cylinder in a heaving and a torsional mode [3]. Through the comparison with the experimental data, the vortex-induced instability of rectangular cylinders was investigated. Moreover we examined the necessity of the three-dimensional simulation of the aeroelastic problem [4]. On the basis of the 0167-6105/95/$09.50 © 1995 ElsevierScienceB.V. All rights reserved SSDI 0167-6105(94)00044-E

T. Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

240

comparison between the three-dimensional computational results of the fixed and the forced-oscillating bluff cylinders, it is confirmed that the flow structure around the cylinder becomes two-dimensionalized due to the oscillation. This means it is expected that the three-dimensional computation is not essentially required for the prediction of the aeroelastic phenomena. Here for the understanding of the fundamental aspects of the aeroelastic instability of a rectangular cylinder, we apply the direct finite difference technique [-5] of the incompressible Navier-Stokes equations to the two-dimensional unsteady flow around a forced-oscillating rectangular cylinder in a torsional mode. The limitation of the two-dimensional simulation is investigated for the prediction of the oscillation problems. With regard to the aerodynamics of an oscillating rectangular cylinder, the relation of the vortex motions and the unsteady aerodynamic forces is examined and we discuss the mechanism of the unstable oscillation and its occurrence by analyzing the pressure distribution on the side surface. The characteristics of the energy distributions are also discussed in order to investigate the vortex-induced and the self-excited instabilities of rectangular cylinders.

2. Problem formulation

The governing equations are given by the following continuity and Navier-Stokes equations, div u = 0 ,

(1)

Ou/Ot + u . grad u = - grad p + 1/Re Au,

(2)

where u, p, t and Re denote the velocity vector, pressure, time and the Reynolds number, respectively, non-dimensionalized by Uo (reference velocity), B (reference length), p (density) and v (kinematic viscosity). For the motion of a body, a moving grid system is incorporated [6]. A time-space transformation of (t, x, y, ..- ) to (~, ~, r/, ... ) is imposed to the original governing equations. The acceleration term of the two-dimensional Navier-Stokes equations is represented as follows, 1 0. + u0~ + vOy = 0~ + j [ ( u -

xdy,

-

1 + j [ - (u - x.)y~ + (v - y . ) x ~ ] 0 , .

(v -

y~)x,]a~

(3)

where J is the Jacobian. This relation means that the convection velocity in the computational domain is replaced by the relative velocity with respect to the moving coordinate. It should be noted that the grid system is moving in the physical domain, but stationary in the computational domain. The governing equations in the

T. Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

241

time-space generalized coordinate system are finally given as follows, x=x(~.q.r),

y=y(~.q.r).

(4) 1

u~ +

[(u -- x ~ ) y , - (v -- y ~ ) x , ] u ¢ + 7 [ -- (u -- x~)y¢ + (v -- y ~ ) x ~ ] u ,

-

j(y.p¢-y,pn)+R~eT(U

,

1

-

l(xep.-x.pe)+R---~7(v,

(5)

7~A = (~A¢~ - 2 f l A ¢ . + 7 A . n ) / J 2 + [(~x~ - 2 f l x c . + y x , . ) ( y e A . - y , A ¢ ) + (~y¢¢ - 2fly¢. + 7Ynn) ( x . A ¢ - x ~ A , ) ] / J 3, J=xcy,--x,y¢.

~=x,

2

+ y , .2

fl=xCx,+ycy,.

(6) 7=x~+Y~.

(7)

/~p = -- { [y,(u¢ -- x¢~) -- y e ( u , -- x,~)] (y, u¢ - y¢ u,) + [xdu.

-

x.d

-

x.(u¢ -

x~)] (y. v~ - y~ v.)

+ (x¢ u, -- x , u ¢ ) [ y , ( v ¢ -- y¢~) -- y¢(v, -- Yn~)]

+ [x¢(v, -- y,~) -- x,(v¢ -- y¢0] (x¢ v, -- x. v¢)}/J 2 + ( y , u¢ - y¢ u, + x¢ v, -- x , v ¢ ) / J A t .

(8)

In order to overcome the numerical instability at high Reynolds numbers, the third-order upwind scheme is employed for the nonlinear convection terms [5]. (U O~x)i ~ Ui -- Ui+ 2 q- 8(Ui+1-- Ui-1) -~" Ui-~- luil Ui + 2 -- 4ui + l + 6ui - 4 u i - 1 + u i - 2 4~x '

(9)

where 6x is the grid spacing. The numerical procedures are based on the MAC method [7]. The pressure field is determined by solving the Poisson equations and the velocity field is computed by the temporal integration of the Navier-Stokes equations.

3. Computational model Figs. la and lb show a computational model where a rectangular cylinder with the depth-breadth ratio ( D / B . D is the depth of the cylinder. B the breadth) equal to 5.0 is

T. Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

242 a

p=0.C u = 1.t3

Fig. 1. Grid system for a rectangular cylinder (400 x 100 points). (a) Whole mesh and the 19ounOary conditions; (b) near the cylinder.

T. Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

243

Fig. 2. Instantaneousvorticitycontours around a stationary rectangularcylinder(D/B = 5.0, Re = 104).

placed in the computational domain (the diameter is 60B). The grid of 400 :x 100 = 40000 points is used for the Reynolds number equal to 104. The approaching flow is assumed to be uniform (U0). For the boundary condition of the pressure around the computational domain, the Dirichlet condition is composed. On the surface of the cylinder, the velocity of the moving grid at the boundary is given and the Neumann condition is employed for the pressure.

4. Torsional oscillation of a rectangular cylinder We present the computational results of the flow around a forced-oscillating rectangular cylinder in a torsional mode (00 = 2 °- i 0 °, where 00 is the amplitude of the torsional oscillation). The displacement angle 0 of the torsional vibration is given as follows, O(t) = Oo sin (2ref,,t),

(lo)

where fm is the frequency of the forced oscillation. Fig. 2 and 3 illustrate the instantaneous vorticity contours around a rectangular cylinder at rest and Vr = 5.0, 10.0 and 30.0 (Uo/f,,B, the reduced velocity). In order to examine the vortex dynamics in the wake of the cylinder, the cases at 2refd = 0 and re (the displacement angle equals 0.0) are shown. The previous experimental data [8, 9] show that the torsional response becomes larger suddenly above Vr --- 5.0 and the resonant phenomenon is observed approximately at Vr = 5.0-8.0. The computational result indicates that the vortex pattern is symmetrically the same and the vortices shed from the frontal comer are merging with the ones shed from the leeward corner at Vr = 5.0. According to the previous experimental results [10], the vortex merging is recognized under the condition that the vortex-induced vibration occurs. The process of vortex merging makes us presume that the vortex excitation appears at Vr = 5.0. In case Vr = 10.0 and 30.0, we can also see the same symmetric patterns of the vortex motions. It is recognized that the cylinder oscillation completely controls vortexshedding whose frequency is coincident with the torsional frequency. According to the experimental data [ 11], the self-excited oscillations based on the aeroelastic instability of the elongated rectangular cylinder (D/B = 4.0) are recognized around V~ = 30. It can be found in the computational result at Vr = 30 that the vortices shed from the frontal

244

T. Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995)239-250

corner never flow leeward constantly, but form the bubble region close to the cylinder surface. These vortices gradually become very strong due to vorticity concentration. As the torsional amplitude becomes larger, the vortices become stronger. The phase-averaged pressure distributions acting on the side surface of an oscillating cylinder at various Vr are shown in Fig. 4. In the case of the response phase (2ztfmt) equal to 0, the effect of the torsional motion on the pressures is discussed in

/d0~ 2rcfmt = O.

2 x f m t = O.

2rcfmt = 1.0~:

2 ~ f m t = 1.0rt

(a) 00= 5 °, Vr=5.0

(b) 00= 5 °, Vr=10.

2rcfmt = O.

2~frnt = O.

2~frnt = 1.0/1;

2rCfmt = 1.Oft

(c) 00= 2 °, Vr=30.

(d) 00= 5 °, Vr=30.

Fig. 3. Instantaneous vorticity contours around an oscillating rectangular cylinder in a torsional mode (D/B = 5.0, Re = 104).

Z Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

-

-

7

245

,

2rcfmt = O.

2rgfmt

1.0rt

=

(e) 00=10 °, Vr=30. Fig. 3. Continued.

cp

Cp

-3

-3

......1 '

(D -2.5 O

~

J! II

Jl--°-'Atrest

-2.5 -2

-1.5

-I .S

~ -0.S ~ 0 ~ O.S

-0.5

:'....I.......i.......L......]L-....I,...,..I .... ........

-I 0 0.5 ....... i . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

-2.5-2-I.5-I

-0.5 0 O.S I

1.5 2 2.5

i ! i ! ~ i ! ! ! ........ ~........ L . . . . - . ; . . - - . . $ . - . . . . ~ . . . - - . . ~ - -........ ...-~.---.~

- 2 . S - 2 - I . 5 - I -0.S 0 0.5 I

1.5 2 2.5

X/B ¢~

.........

i}iii!!i!

I

X/B

Co

cp

1

1

o.s 0

t ~ -0.S

l i I ......

÷

.............

i- -.0"-

.~--

-. • ---~ .-------

.-

-'-'i

........

i ........

" .......

0

i

i

i

-O.S

i

i

i q ~ ' i - - ~ a - ~

~

r~

~ -1 .S

-1.5

N -z

-z

© -2.5

.......]-.

-2.s [ I

-3

-2.S-Z-l.S-I

-0.5 0 0.5

1 l.S

2 2.5

W=~o.o,~m~oS.o(ao~

. . . . . A, r~,

-3 -2.5-2-I.5-I

-0.5 0 O.S 1 1.5 2 2.5

X/B

X/B

(a) 00=5 °, Vr=5.0,

...... L.......

(b) 0 0 = 5 °, Vr=lO.

Fig. 4. P h a s e - a v e r a g e d pressure distribution acting on the side surface of an oscillating rectangular cylinder at 2nfmt = 0.0. in a torsional m o d e (D/B = 5.0, Re = 104).

co -3

- -

-Z3

co

• Vr-30.0,Amp-Z.0deg --A--Exp. by Shiraishi. et el.

J....

(Vr-Z4.Z,Amp-4.46deg) |

!ii

,~

7jk:°::~L~'! ....................................... !~

-O.S 0 0.5 1 1.5 2 2.5

-1 .S -1 -0.S 0 i i ...............i........L i i _ l 03 I -2.5-2__-1.5-1 -0.5 0 0.5 1 1.5 2 2.5

X/B

X/B

-1.5 ~

• Vr-30.0,Amp-S.0deg L -3 m r [ . .......[J --A-- Exp. by Shiraishi. et el. L.J -2.5 '-'"'il (Vr-Z4.Z,Amp-4.46deg) 1...I -Z

~--

-0.5 i i .."~] ' ' ~ -( ~ i i i i i i i i ....... ~.......i..-..-..ii........~....---*.--...-i........il . . . ........ ....... 0.5 1

-2.5-2-1.5-1

l

co

co 1

¢~

~

0.5 -0.5 -1 -1.5

1

;
-1.5

i

i

i

i

i

i

- - + - = - •

.------}--[

a"

'

i

i

i

•

Vr-30.0,Amp-S.0deg by Shiraishi. et el. (Vr-Z4.Z,Amp-4.46deg) - -o-- At rest

-2 ~--A"Exp.

-2.5

(Vr=Z4.Z,Amp-4.46deg) - - o - - A t rest -0.S 0 0.S i

i

~

-o.s ~

iiiiiiiii:i iiiiill " 0.0,Am0-2.0,e0 Vr........................... --A--Exp. by ShiraishL et el. -2.5-2-1.5-1

i

o

.......
~-2.5 ~ -3

i

-3 -2.5-2-1.5-1-03

1.5 2 2.S

0 0.5 1 1.5 2 2.5

X/B

X/B

(c)~Oo--2 o, Vr=30.

(d) 00=5 °, Vr=-30.

Cp "~ •

Vr-30.0,Amp-10.0deg h -3 ---.-.-I--A--Exp. by Shiraishi. et el. I.----I ~.~-2.5 ---.-.] (Vr=Z4.2,Amp=4.46deg)[.....i -1.5 ................,.........................

=

~

~. -O 5~ i

'

i

]

I

i

~i

i........i

.

i

i _ ! ........

i........L.......L......i_..] ......

ili

o.5

.

!ii

...........................................................................

1 -2.5-2-1.5-1

-0.5 0 0.S 1 1.5 2 2.5

X/B [- =

Vr=30.0,Amp=10.0deg by ShiraishJ. et el. 1 i] (Vr=24.2,Amp=4.46deg) 0.5 ........iL . . . . . At rest

Cp.~j --A--Exp.

O tJ

0 -0.5

i

i

~

- 1 ~'~'~ "i5"~"~ -1 .S

© ~

-2 -2.5

O

-2.5-2-1.5-1

i •

:

io-

=

i

.......

~

i_..L.

~'"""i

i

-0.S 0 0.S 1 1.5 2 2.5

X/B

(e) O0=lO °, Vr=30. Fig. 4. C o n t i n u e d .

T. Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

247

comparison with the case of a stationary cylinder. It is recognized that when V~ = 5.0 and 10.0 the pressure on the lower surface becomes negative after separation, changing into positive in the leeward region. On the other hand, when Vr = 30 the negative pressure recovers slowly and still remains in the same leeward location, where the self-excited force is induced. In detail, we can see a few peaks of negative pressure, which cannot be seen in the previous experimental data 1-12]. This discrepancy must be due to the two-dimensional computation, however this results never have a significant effect on the evaluation of the aerodynamic characteristics. Approximately it can be presumed that the present simulations show qualitatively good results, for example the peak location of the negative pressure on the lower surface is moving leeward, as V~ becomes higher. The results for the various amplitudes of 2 °, 5 ° and 10 ° at V~ = 30 present almost the same tendency. As the amplitude becomes larger, the pressure fluctuations also become larger. The unsteady characteristics of the aerodynamic forces are investigated by analyzing the energy distribution on the side surface of the cylinder in Fig. 5. The total energy distribution is given as follows, Etotal =

Cp dA r dO,

(11)

where Tp, A, Cp, r are the reduced natural period, area, pressure coefficients, arm of the moment caused by each Cp. The energy distributions indicate the consistent result with the pressure distribution on the side surface in Fig. 4. With regard to the effects of Vr, it is recognized that the self-excited forces mainly appear in the windward region at Vr = 5.0 and 10.0, and in the leeward region at V, = 30. The aerodynamic damping force is very large in the leeward region when at Vr = 5.0 and 10.0, so finally the total torque on the rectangular cylinder becomes a damping for the oscillation. With regard to the dependence of the torsional amplitude at V~ -- 30, the results have the same tendency with two peaks appearing on the side, one immediately after the frontal corner and the other downstream of the center of torsion. Due to these two peaks, the cases at V~ = 30 become aeroelastically unstable conditions like the torsional flutter. Through the comparison with experimental data [12], the computational results are in good agreement with the experimental data concerning the energy distributions based on the pressure as well as the aeroelastically unstable conditions of the oscillating cylinder. Time variations of the torque coefficients (Cu) are given for several displacement angles of an oscillating cylinder in Fig. 6. In case of Vr = 5.0, the value of CM is smoothly changing as a sinusoidal curve under vortex excitation, however at V~ = 30.0 the CM curve becomes complicated because the vortices are evolving and changing their shapes close to the side surface. The phase angle of the torque to the oscillation becomes negative at Vr ----5.0 and 10.0, it can be confirmed that this case aeroelastically becomes stable. On the other hand, the phase becomes positive at Vr = 30.0 and we can also confirm the energy transfer of the divergent oscillation. This tendency is in good agreement with the experimental data [11].

T. Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

248

>,0.3

" ~.

'

o

.

1

~

" ~

o

°~-2.5-2-1.5-1,0.5 0 0.5 1 1.5 2 2.5 x/e

(a) Computationalresults(00=5°.) >,

~

F~----~-- Amplitude = 2.0 d e g r e e Amplitude = 5.0 degree m--- Amplitude = 7.0degree

0.15

.

.

.

.

0.1

e-

LU 0.05 t1:l

-0.05 -0.1 -0.15 - 2 . 5 - 2 - 1 . 5 - 1 - 0 . 5 0 0.5 1 1.5 2 2.5 X/B

(b) Computationalresults(Vr=30.) 0.15 fLU

0.1 0.05

I

I

o

I

I

I

L

I

I

exl)

/

i

i

,.,i ......... i.........i......... - .........i......... i.........i .........i......... i • -4-.%

l-

I

----e-.- Exp(Vr = 24.2, Amp = 4.46 degree) I

i

i

~

4

:

0 ca) -o.o5 ..................i......................................................................i -0.1 ................. i............................................................................. .

-0.15 -2.5-2-1.5-1-0.50

.

.

.

.

.

.

.

i 0.5 1 1.5 2 2.5 X/B

(c) Experimentalresults by Shiraishi et al. ( 00=4.46

°,

Vr=24.2.)

Fig. 5. Total energy distribution on the side surface of an oscillating rectangular cylinder in a torsional m o d e (D/B = 5.0, Re = 104).

T. Tamura et al./J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

81

CM I

4

i

[ . . . . . ,o,l,a,) I

I

.,-......, ........ , . . . .

I

i

I

,.......4 ........ ,...

0.05 0 o. ;T. -0.05 -0.1 -0.15 -0.2 176.25

~0

166.25

o.1

, ........ ~ ......

: ~ 2

-8 161.25

0.2 0.15

I

171.25

tUO/B

(a) Vr=5. 8 6 4 2

-2 -4 -6

~[ CM][ I ..... r°'
g..,

.................

i--,

........................................

i ..........

i

i

70 170

60

180~

0.2 0.15 0.1 0.05 0 -0.05 -0.1 -0.15

19d 1 °2

tUOIB

(b) Vr= 10. 8

6 J-f ".

4

c.

_

_ _

,"~ ..........

! ..........

t [ .....

0.2 0.15

rot(rad)

'T"".

........ , .......... ; .......... , .......... , .......... 4 .......... , .......... $ .......... i ......... !

i

~0

.!

i "~i

"2

-6

Ir

~

......... i ..........

-8

400.00

i

!

m.,:

!a

,~ "

i

i

i

i

i

•

~

ikr

i

i

"1/

.

o.1 0.05 ~"

i

,

o

-0.05 -0.1 ~ .......... ; .......... i .......... ~ .......... i .......... ~..................... -0.15 i i J i i i -0.2 430.00 460.00 490.00 !

m:

!

~

m]

(c) Vr=30.

i

tUO/B

Fig. 6. Time variations of the torque coefficients (CM) and the torsional angle (0o = 5°).

249

250

T. Tamura et aL /J. Wind Eng. Ind. Aerodyn. 54/55 (1995) 239-250

5. Conclusion W e s i m u l a t e the u n s t e a d y flow a r o u n d a f o r c e d - t o r s i o n a l l y - o s c i l l a t i n g r e c t a n g u l a r c y l i n d e r b y m e a n s of the direct finite difference t e c h n i q u e w i t h o u t a n y t u r b u l e n c e model. T h e flow p a t t e r n s a r o u n d the c y l i n d e r s h o w the distinctive characteristics such as a v o r t e x m e r g i n g in v o r t e x e x c i t a t i o n o r a v o r t e x c o n c e n t r a t i o n in t o r s i o n a l flutter. T h e u n s t e a d y forces with r e g a r d to the a e r o d y n a m i c instability are discussed focusing o n the energy d i s t r i b u t i o n s on the cylinders. A c c o r d i n g l y it is s h o w n t h a t the p h y s i c a l m e c h a n i s m of u n s t a b l e oscillations like a v o r t e x - i n d u c e d v i b r a t i o n a n d a t o r s i o n a l flutter are clarified b y a n a l y z i n g the n u m e r i c a l d a t a with respect to the pressure distributions.

References [1] T. Tamura, E. Krause, S. Shirayama, K. Ishii and K. Kuwahara, Three-dimensional computation of unsteady flows around a square cylinder, in: Proc. 1lth ICNMFD, Williamsburg, 1988. 1-2] T. Tamura and K. Kuwahara, Numerical analysis on aerodynamic characteristics of an inclined square cylinder, AIAA Paper, 89-1805, 1989. [3] T. Tamura and K. Kuwahara, Numerical study on aerodynamic instability of oscillating rectangular cylinders, J. Wind Eng. Ind. Aerodyn. 41 (1992) 253, 254. I-4] T. Tamura, Y. Itoh, A. Wada and K. Kuwahara, Numerical investigation on the aerodynamic instability of bluff cylinders, 1st CWE 92, J. Wind Eng. 52 (1992) 557-566. [5] T. Tamura, K. Tsuboi and K. Kuwahara, Numerical simulation of unsteady flow patterns around a vibrating circular cylinder, AIAA Paper, 88-0128, 1988. [6] T. Kawamura and K. Kuwahara, Computation of high Reynolds number flow around a circular cylinder with surface roughness, AIAA paper, 84-0340, 1984. [7] F.H. Harlow and J.E. Welch, Numerical calculation of time-dependent viscous incompressible flow of fluid with free surface, Phys. Fluids 8 (1965) 2182-2189. [8] Y. Nakamura, M. Nakashima and K. Watanabe, Bulletin R.I.A.M., Kyushuu Univ. 60 (1984) 45-52 [in Japanese]. [9] T. Yoshimura, A. Miyake, T. Hirayama and T. Akamatsu, Vortex excitation of elongated rectangular prism, in: Proc. Symp. on Wind engineering, 1986 [in Japanese]. [10] N. Shiraishi and M. Matsumoto, On classification of vortex-induced oscillation and its application for bridge structures, J. Wind Eng. Ind. Aerodyn. 14 (1983) 419430. [11] H. Tomizawa, K. Washizu, A. Ohya and Y. Otsuki, Wind tunnel experiments of box-like models of buildings: Measurement of unsteady aerodynamic forces by forced oscillationmethod, Reports of the Research Laboratory of Shimizu Corporation, No. 20, 1973 [in Japanese]. [12] N. Shiraishi, M. Matsumoto, H. Shirato and T. Memita, On unsteady pressure characteristics of rectangular cylinders in flutter, in: Proc. 9th Nat. Symp. on Wind engineering, 1986 [in Japanese].