Numerical Simulation of the Effect of a Single Gust on the Flow Past a Square Cylinder

Abstract

The flow past a square cylinder under the influence of a one dimensional gust was investigated using computational fluid dynamics (CFD). The effect of upstream wind gusts of the same amplitude but different duration was investigated with respect to their effect on the flow, the vortex-shedding, and the pressure distribution around the square cylinder. For the computations, a very large eddy simulation (VLES) model was implemented in an in-house code and validated against numerical and experimental results from the literature. The gusts of different duration were found to have a distinctly different effect. The short-duration gust causes a lock-on behavior with cessation of the alternating vortex shedding, and a symmetric pair-vortex was created above and below the square cylinder. It was observed that the pressure distribution on the lateral sides of the cylinder has the same magnitude and phase, which resulted in a zero total lift coefficient. In terms of a free-standing structures, such as a building, this would lead to zero instantaneous forces and pressure difference in the lateral direction with obvious implications for dynamic response and cross ventilation.

1. Introduction

Turbulent and non-stationary winds are often observed in the atmosphere, having the common characteristic of different speeds at different timescales [1]. Large fluctuations in wind speed over small timescales are associated with wind gusts [2] and are of interest to several applications such as wind turbine operation, loads on structures, aeronautics [3,4], and rail transport [5].

Most studies of wind-induced gusts have been experimental ones that measured the flow characteristics in, e.g., urban environments throughout the world, such as Newcastle, Australia [6], or Shangai and Shenzen, China [7,8]. Focus is usually placed on the yet-to-be well-understood urban flow physics [7,8] but also on the characterization of the gust itself [6]. Even outside urban areas, research into the meteorological physics of the gusts is active for local phenomena such as the Vento Norte in Southern Brazil [9] or coastal areas, e.g., Lake Victoria [10] or Pearl River Estuary [11]. Yet another active area of research is extreme gust phenomena such as typhoons [12,13] or tropical cyclones [14]. Numerical studies are not as widespread, and most concern wind turbine operation during wind gusts. Researchers have found ways to reduce loads exerted on wind turbines [15,16] and to make them as efficient as possible when they are subjected to extreme gusts by designing and tuning fuzzy logic controllers [17] and optimizing the design of the wind turbine blades [18,19].

On the other hand, only a limited number of numerical studies of an unsteady, gusty flow past a bluff body have been published. Some examples are a sinusoidal signal past a circular cylinder [20] or a single coherent gust past a flat plate [21], while other studies have used gust inlet velocities that describe a sequence of gusts [22], possibly even of a stochastic nature [23,24], and are more similar to oscillation functions than to a single gust described by U.S. weather-observing practice [25].The most common design standards prescribed to provide guidelines in order to take wind gusts into account and to satisfy U.S. weather-observing practice [25] are Kasperski’s standard [26] and the standard defined by the International Electrotechnical Commission (IEC) 61400.2-2013 [27]. The latter describes the extreme operating gust (EOG), which is a representative gust profile forming a symmetric Mexican hat shape. This standard is used to describe the gust for a wind turbine during its operation, and this is probably the reason that it has not been used to describe gusts to which bluff bodies are subjected [28].

Turbulent flow past bluff bodies has attracted significant attention in recent years, with a cylinder of square cross section being one of the most common geometries studied, where vortex shedding, flow separation, and shear-layer instability appear in the vicinity of the bluff body [29]. The flow past a square cylinder may be considered relevant to many engineering applications such as modern buildings, towers, high-rise skyscrapers, and vertical columns of long-span bridges because they share similar aerodynamic characteristics [30]. Many studies have been carried out for this type of flow, investigating wake vortex dynamics [31], the interaction with the ground [32], dependence on Reynolds number [33], the effects of 3D flow structures [34], a universal Strouhal number [35], and vortex-shedding suppression [36].

The complex flow separation, vortex-shedding, and recirculation zones on and in the wake of the square cylinder may induce large instantaneous aerodynamic forces, and the investigation of the drag and lift forces exerted on its external surfaces provides important information for numerous applications. One of the less investigated applications in this context is that of the operation and energy performance of a building and of occupant comfort (natural ventilation, infiltration, wind loads), in spite of its obvious dependence on the pressures exerted on the outer surfaces. It is well-established that infiltration and ventilation is a crucial factor when it comes to sustainable energy buildings [37] and healthy environments [38], and many relevant studies have been carried out regarding the impact of opening positions [39], the effect of neighboring objects [40], and also the accuracy of applied numerical methodologies [41].

However, to the authors’ knowledge, the effect of single, instantaneous fluctuations such as gusts in upstream conditions has not been extensively investigated until now. Thus, the novel aim of the present numerical study is to introduce a realistic wind gust, which satisfies the IEC Standard [27] and U.S. weather-observing practice [25], upstream of the flow past a square cylinder and to compare the results for gusts of varying duration.

Although there are many experimental studies of the flow past a square cylinder conducted, e.g., by measuring static pressure distributions at midspan in a wind tunnel [42] or by directly measuring aerodynamic coefficients in wind tunnel [43] and water channel tests [44,45], the numerical simulation of this flow is challenging. Several numerical studies have been performed using and comparing RANS, URANS, and LES turbulence models [46,47,48,49], with attempts to perform DNS being much less common [50] even if considerably old [51]. Recently, hybrid RANS/LES methodologies have appeared, i.e., DES [52], VLES [53], and PANS [54], which are significantly more attractive in terms of demands on computational resources compared to DNS and LES and may provide a better representation of unsteady turbulent flows than URANS. For the present study, the Han and Krajnovic [55] VLES model is implemented in an in-house code. The implementation is validated against the measurements of Lyn et al. [45] for the steady flow past a square cylinder at Re = 22,000 in order to reliably apply it to studying the effects of upstream gusts. After validation of the code, single gust events are applied to the inlet boundary of this flow, and the effects of two cases of different gust duration are compared. It is found that although the two gust events have the same velocity structure with a different duration, and they both induce high instantaneous pressures on the cylinder’s outer surfaces, their effect on the global indices such as lift, drag, or pressure differences are markedly different. Implications to infiltration and ventilation of buildings are discussed.

2. Materials and Methods

2.1. Numerical Code

The methodology applied in the in-house code for the numerical computation of the flow field solves the Navier–Stokes equations using the finite-volume method. The generalized form of the governing equations is

$$\frac{\partial }{\partial \mathsf{t}}(\mathsf{\rho }\Phi )+\frac{\partial }{\partial \mathsf{x}}(\mathsf{\rho }\mathsf{u}\Phi )+\frac{\partial }{\partial \mathsf{y}}(\mathsf{\rho }\mathsf{v}\Phi )+\frac{\partial }{\partial \mathsf{z}}(\mathsf{\rho }\mathsf{w}\Phi )-\frac{\partial }{\partial \mathsf{x}}({\Gamma }_{\Phi }\frac{\partial \Phi }{\partial \mathsf{x}})-\frac{\partial }{\partial \mathsf{y}}({\Gamma }_{\Phi }\frac{\partial \Phi }{\partial \mathsf{y}})-\frac{\partial }{\partial \mathsf{z}}({\Gamma }_{\Phi }\frac{\partial \Phi }{\partial \mathsf{z}})={\mathsf{S}}_{\Phi }$$

(1)

where u, v, and w are the corresponding velocity components in local Cartesian coordinates, and ρ is the density of the fluid. Φ is a generalized variable that takes values Φ = 1, u, v, and w for the continuity and momentum equations, respectively, or Φ = k,ε for the conservation equation of turbulent kinetic energy (k) or its dissipation rate (ε), respectively. For the continuity equation, S_Φ in (1) is zero, while for the momentum equations, it is written in tensor form as:

$${\mathsf{S}}_{\mathsf{\Phi }={\mathsf{u}}_{\mathsf{i}}}=-\frac{\partial \mathsf{{\rm P}}}{\partial {\mathsf{x}}_{\mathsf{i}}}+\frac{\partial }{\partial {\mathsf{x}}_{\mathsf{j}}}\left(\mathsf{\mu }\cdot \frac{\partial {\mathsf{u}}_{\mathsf{j}}}{\partial {\mathsf{x}}_{\mathsf{i}}}\right)$$

(2)

where x_i is equal to x, y, z and u_i to u, v, w when i is 1, 2, 3 respectively. In addition, for the k and ε conservation equations, the source term S_Φ is:

$${\mathsf{S}}_{\Phi =\mathsf{k}}=\mathsf{\mu }\frac{\partial {\mathsf{u}}_{\mathsf{i}}}{\partial {\mathsf{x}}_{\mathsf{j}}}\left(\frac{\partial {\mathsf{u}}_{\mathsf{i}}}{\partial {\mathsf{x}}_{\mathsf{j}}}+\frac{\partial {\mathsf{u}}_{\mathsf{j}}}{\partial {\mathsf{x}}_{\mathsf{i}}}\right)-\mathsf{\rho }\cdot \mathsf{\epsilon }$$

(3)

$${\mathsf{S}}_{\Phi =\mathsf{\epsilon }}={\mathsf{c}}_1\frac{\mathsf{\epsilon }}{\mathsf{k}}\frac{\partial {\mathsf{u}}_{\mathsf{i}}}{\partial {\mathsf{x}}_{\mathsf{j}}}\left(\frac{\partial {\mathsf{u}}_{\mathsf{i}}}{\partial {\mathsf{x}}_{\mathsf{j}}}+\frac{\partial {\mathsf{u}}_{\mathsf{j}}}{\partial {\mathsf{x}}_{\mathsf{i}}}\right)-{\mathsf{c}}_2\mathsf{\rho }\frac{{\mathsf{\epsilon }}^2}{\mathsf{k}}$$

(4)

In Equations (1)–(4), the other variables appearing are Γ_Φ = μ/σ_Φ, μ = μ_l + μ_t, where μ_l is the dynamic viscosity of air, and μ_t is the turbulent viscosity, with μ_t = C_μk/ε for the standard k-ε turbulence model. In addition, constants are c₁ = 1.44, c₂ = 1.92, c_μ = 0.09, σ_ui = 1, σ_k = 0.9, and σ_ε = 1.3.

The code is a standard finite volume implementation of the SIMPLE algorithm [56] on a Cartesian grid with collocated variables and the discretized equations are solved using a tri-diagonal matrix algorithm (TDMA). For the convective terms, the bounded second-order upwind scheme (BSOU) [57] is used, which combines and bounds the first- and second-order upwind schemes, providing better accuracy than the hybrid scheme and reducing the numerical diffusion problems. Previous application of this CFD model to steady state external flows past buildings and urban geometries has proven reliable [58,59,60].

For turbulence closure, we apply the Han and Krajnovic [55] VLES concept with the k-ε model, where the resolution control function, Fr, was bounded between zero and one:

$${\mathsf{F}}_r=\min \left(1.0,{\left[\frac{1.0-\exp (-\frac{\mathsf{\beta }{\mathsf{L}}_{\mathsf{c}}}{{\mathsf{L}}_{\mathsf{k}}})}{1.0-\exp (-\frac{\mathsf{\beta }{\mathsf{L}}_{\mathsf{i}}}{{\mathsf{L}}_{\mathsf{k}}})}\right]}^{\mathsf{n}}\right)$$

(5)

L_c, L_i, and L_K are the turbulent cut off length scale, integral length scale, and Kolmogorov length scale, respectively, defined as:

$${\mathsf{L}}_{\mathsf{i}}={\mathsf{k}}^{1.5}/\mathsf{\epsilon }\hspace{1em}{\mathsf{L}}_{\mathsf{c}}={\mathsf{C}}_{\mathsf{x}}{(\Delta x\cdot \Delta y\cdot \Delta z)}^{\frac{1}{3}}\hspace{1em}{\mathsf{L}}_{\mathsf{{\rm K}}}={\left(v^3/\mathsf{\epsilon }\right)}^{1/4}$$

(6)

In Equation (6), the model constants are C_x = 0.61 and n = 2, and the recommended value of β is β = 2 × 10⁻³ based on the numerical studies in the original Speziale VLES model [53]. It should be noted that near the wall, L_c > L_i, leading to Fr→1 (see Equation (5)), so the hybrid model recovers to the RANS model, a behavior similar to the (detached eddy simulation) DES concept [52].

Turbulent viscosity from the k-ε model is calculated as:

$${\mathsf{\mu }}_{\mathsf{t}}^{\mathsf{s}\mathsf{u}\mathsf{b}}={\mathsf{F}}_r\mathsf{\rho }{C}_{\mathsf{\mu }}\frac{{\mathsf{k}}^2}{\mathsf{\epsilon }}$$

(7)

2.2. The Flow Field Past a Square Cylinder

Validation is performed for the flow past a square cylinder according to the experiment of Lyn et al. [45]. At Re = 22,000, with air at atmospheric pressure and standard temperature, the cylinder edge length was chosen to be D = 1 m on the XY plane (Figure 1a), the height of the square cylinder is parallel to the z-axis, and the direction of the flow is along the x-axis with an inlet velocity of U_inlet = 0.33 m·s⁻¹.The computational domain is 20D × 14D × 4D, and the center of mass of the square cylinder is at (x,y) = (0,0) (Figure 1a). No-slip wall boundaries simulate the water tunnel in the lateral direction, and free-slip boundaries were assumed for the top and bottom of the domain in the spanwise direction. Zero gradient boundary conditions were imposed at the outflow, which was 15D away from the rod center, and no pressure reflections were observed. In this paper, the flow was examined with a mesh (Figure 1b) that consists of 0.378 million cells (M₃₇₈) in total. This is similar to the coarsest (M₁) mesh that Han and Krajnovic used and found to be sufficient. The grid is clustered near the wall, and the first node is located at Δy/D = 2.5 × 10⁻³. In the present study, grid independency was also investigated with a finer mesh (M₅₇₅), which consisted of 0.575 million cells. Further details on the performance of the numerical setup are presented in Section 3.1.

The non-dimensional time step Δt* = Δt/(D/U_inlet) is equal to 0.11. At the inlet, the mean horizontal velocity U is uniform and follows the horizontal direction, the turbulent intensity is set to 2%, and the turbulent viscosity ratio is equal to 20.0.

2.3. Gust Introduction

According to U.S. weather-observing practice [25], a gust is a sudden, brief increase in the speed of the wind, which lasts less than 20 s. The peak of the wind speed should reach at least 16 knots, and the variation in wind speed between the peaks and lulls should be at least 9 knots [25]. Although gusts of shorter duration (~3–5 s) are most common, gusty turbulent fluctuations may also last in the order of minutes [61] also in complex terrain [62].

For numerical simulation, some form of gust function must be applied to the inlet velocity profile of the above-mentioned flow field. Standardization of gust definition has been performed by IEC work groups, and one of the extreme wind conditions, described by the International Electrotechnical Commission (IEC) standard IEC 61400.2-2013 [27], is the extreme operating gust (EOG). This representative gust profile agrees with gust events that were found to occur 100 times during a 12 month interval [6,27] and takes on the form of a symmetric Mexican hat shape with a decrease in speed followed by a steep rise, a steep drop, and a rise back to the original value. The IEC standard is usually applied to describe the gust for a wind turbine in operation, but in the present study, it was also considered valid in order to define the gust velocity profile upstream of the bluff body being studied:

$$\mathsf{u}(\mathsf{t})=\left\{\begin{array}{c}\mathsf{u}(\mathsf{z})\\ \mathsf{u}(\mathsf{z})-0.37\cdot {\mathsf{V}}_{\mathsf{g}\mathsf{u}\mathsf{s}\mathsf{t}\mathsf{N}}\cdot \sin \left(\frac{3\pi \mathsf{t}}{\mathsf{T}}\right)\cdot (1-\cos \left(\frac{2\pi \mathsf{t}}{\mathsf{T}}\right)\\ \mathsf{u}(\mathsf{z})\end{array}\begin{array}{c}\mathsf{t}<0\\ 0\le \mathsf{t}\le \mathsf{T}\\ \mathsf{t}>\mathsf{T}\end{array}\right\}$$

(8)

In Equation (8), a normal steady-state profile u(z) is assumed to persist up to t = 0 s when the gust begins. The gust period is the total event duration, considered to be T = 10.5 s for a recurrence period of N = 1, and V_gustN is the gust amplitude at the hub height, which is the height of interest in the case of a wind turbine [6]:

$${\mathsf{V}}_{\mathsf{g}\mathsf{u}\mathsf{s}\mathsf{t}\mathsf{N}}=\mathsf{\beta }\cdot \frac{{\sigma }_1}{1+0.1\left(\raisebox{1ex}{$\mathsf{D}$}\!\left/ \!\raisebox{-1ex}{${\mathsf{\Lambda }}_1$}\right.\right)}\hspace{1em}{\mathsf{\Lambda }}_1=\left\{\begin{array}{c}0.7\cdot {\mathsf{z}}_{\mathsf{h}\mathsf{u}\mathsf{b}}\\ 21 \mathsf{m}\end{array}\begin{array}{c}\mathsf{w}\mathsf{h}\mathsf{e}\mathsf{n}\\ \mathsf{w}\mathsf{h}\mathsf{e}\mathsf{n}\end{array}\begin{array}{c}{\mathsf{z}}_{\mathsf{h}\mathsf{u}\mathsf{b}}<30 \mathsf{m}\\ {\mathsf{z}}_{\mathsf{h}\mathsf{u}\mathsf{b}}\ge 30 \mathsf{m}\end{array}\right\}\hspace{1em}{\sigma }_1={\mathsf{{\rm I}}}_{15}\frac{15+\alpha \cdot {\mathsf{V}}_{\mathsf{h}\mathsf{u}\mathsf{b}}}{\alpha +1}$$

(9)

where Λ₁ is the turbulence scale parameter, β = 4.8 is the parameter for the EOG model when N = 1, D = 5 m is the rotor diameter, and σ₁ is the standard deviation. I₁₅ = 0.18 is the dimensionless characteristic value of the turbulence intensity at 15 m/s, and α = 2 is the dimensionless slope parameter for a class III SWT. Here, V_hub is considered to be the inlet mean horizontal velocity.

The gust function and parameters in Equation (8) refer to a full-scale velocity profile and gust period (T_n). For the present situation, this may be reasonably considered to correspond, e.g., to a full-sized square cylinder (e.g., building structure) of side L_n = 10 m and an upstream wind velocity of u_n = 6.3 m/s. In order to facilitate comparison with the results of the steady upstream velocity validation test presented above, the gust function should be scaled and a corresponding gust period (T_m) be calculated. Based on standard scaling procedures [63], we can define a time-scale similarity ratio between fullscale and validation (model) scale:

$$\frac{{\mathsf{T}}_{\mathsf{m}}}{{\mathsf{T}}_{\mathsf{n}}}=\frac{\raisebox{1ex}{${\mathsf{L}}_{\mathsf{m}}$}\!\left/ \!\raisebox{-1ex}{${\mathsf{u}}_{\mathsf{m}}$}\right.}{\raisebox{1ex}{${\mathsf{L}}_{\mathsf{n}}$}\!\left/ \!\raisebox{-1ex}{${\mathsf{u}}_{\mathsf{n}}$}\right.}=1.9$$

(10)

In the validation test case, the square cylinder had a side length of L_m = 1 m, and the inlet streamwise velocity was u_m = 0.33 m/s. Thus, if the full-scale gust period is T_n = 10.5 s (in agreement with the standard of EOG [27]), then the gust period of the model square cylinder would be T_m ≈ 20 s. From Equations (8) and (9), the model gust velocity signal is calculated and presented in Figure 2. The same Figure also shows a wind gust of the same amplitude but an arbitrarily exaggerated model scale duration of 100 s, corresponding to a full-scale gust period (Equation (10)) of 52 s (a little less than one minute), which was also examined for comparison purposes.

3. Results and Discussion

Validation of the VLES implementation in the code was performed by comparison to results from the literature, and then, the gust events described previously were introduced to the upstream flow and the flow field around the cylinder, and the pressure distribution on its outer surfaces were calculated.

3.1. VLES Validation of Flow Past a Square Cylinder

Validation of the VLES implementation in the in-house code was performed by comparison to the experimental measurements of Lyn et al. [45] and the numerical calculation of Han and Krajnovic [55] for the turbulent flow past a square cylinder. The computational domain and the boundary conditions have already been discussed in Section 2.2. However, for the study of the wind gusts, it was also necessary to examine this flow when free-stream boundary conditions are imposed in the lateral direction instead of the wall boundaries present in the water channel realization.

Calculation results of some of the flow field global parameters are presented in

Table 1. Comparison of the global flow parameters including results from Sohankar et al.’s LES [47], Han et al.’s VLES [55], Lyn et al.’s experiment [45], and Baroneet al.’s DES [52].

	Re	<CD>	CDrms	CLrms
Present VLES-M378 (No-slip wall)	22,000	2.28	0.22	1.62
Present VLES-M378 (Free-slip wall)	22,000	2.25	0.27	1.61
Present VLES-M575 (No-slip wall)	22,000	2.33	0.22	1.64
VLES -M1 [55]	22,000	2.21	0.18	1.36
VLES -M2 [55]	22,000	2.28	0.19	1.51
LES [47]	22,000	2.03–2.32	0.16–0.20	1.23–1.54
DES [52]	19,400	2.11	0.26	1.16
Experiment [45]	21,400	2.10	-	-

and compared with previous studies, and grid independency is also verified. Similar mean and RMS drag and lift coefficients are computed with theM₃₇₈ and M₅₇₅ meshes, giving a maximum deviation of ~2% for <C_D>.

Minimal deviations are observed when free-slip and no-slip boundaries were set, and overall, the results are in agreement with those of other studies given the rather large spread in those results. An exception is the standard deviation of the lift coefficient, which is close only to the lift coefficient computed with the denser grid by Han and Krajnovic VLES [55]. It has to be mentioned that the precise grid distribution of the meshes in other studies is not known. Given the level of agreement and the spread in the results among other studies, it has been considered that mesh M₃₇₈ is adequate, and all the calculations below are performed on this mesh.

In Figure 3, the timeseries of the streamwise velocity component at the centerline, five cylinder side lengths downstream (x/D = 5), is shown for 234 non-dimensional time units, which correspond to ~30 complete natural shedding periods. The frequency of these periods gives a Strouhal number equal to St = f.D/U_inlet = 0.122, which is the same as the value calculated in [55].

In terms of the behavior of the turbulence modeling, the mean y⁺ values of the first node along the square cylinder sides are presented in Figure 4 and shown to fall within the range of y⁺ = 1 − 4.5. Figure 5 shows the ratio of the effective to laminar viscosity as an indication of the performance of the VLES model. Effective viscosity is limited to an order of magnitude increase compared to laminar viscosity, mainly in areas of high shear (shear layers from cylinder sides and vortex cores), but for the most part of the flow, the increase is about fivefold compared to orders of magnitude commonly found in URANS studies.

The mean flow field is computed for the above-mentioned duration of 234 non-dimensional time units and the mean streamwise velocity profile along the centerline behind the cylinder is presented in Figure 6. Results of the present study are in good agreement with Lyn et al.’s experiment [45] and also agree with those of Han and Krajnovic [55], with a small divergence of the order of 4% far downstream, where RANS behavior is recovered, and grid spatial resolution is coarser. Franke and Rodi’s [64] calculation using the standard k-ε model shows an unrealistically large recirculation zone, while Sohankar and Davidson’s LES [47], also shown for comparison, overestimated flow recovery in the wake.

The profiles of the three fluctuating (RMS) resolved velocity components along the central line of the wake are shown in Figure 7. The u_RMS velocity is in good agreement with both the experimental measurements and Han and Krajnovic’s VLES [55] within the recirculation zone behind the cylinder but fails to predict the behavior in the far wake. However, the v_RMS velocity component is in excellent agreement with the experimental measurements, while it was underestimated by Han and Krajnovic’s VLES [55] and LES [47]. This should also be considered in relation to the higher RMS value of the lift coefficient in Table 1. Experimental data were not available for the w_RMS velocity component, but the predictions of the models are generally in close agreement to each other. Total fluctuation energy (TFE) calculated along the central line of the wake is shown in Figure 8. The total fluctuation energy is computed from the three fluctuating resolved velocity components as TFE = $(\overline{u{\prime }^2}+\overline{v{\prime }^2}+\overline{w{\prime }^2}$)/2 and is compared to Lyn et al.’s experiment [45] and to Bouris and Bergeles’s 2D LES results [65]. Finally, present VLES is in good agreement with the experimental measurements within the recirculation zone and the near wake, but the effect of the u_RMS velocity deviations observed in Figure 7 becomes evident in TFE in the far wake, i.e., downstream of x = 3D.

Mean velocity profiles in the transverse direction at x/D = 1 are presented in Figure 9. It is verified that in the near wake, both mean streamwise and transverse velocity profiles are in good agreement with the experiment.

Furthermore, streamlines of the average flow are presented in Figure 10 and agree well with the results of Trias et al. [50], who conducted a DNS study past a square cylinder at Re = 22,000. The two large recirculation zones appearing at the top and leeward sides of the cylinder (Figure 10) are of the same overall size as those predicted by the DNS method [50] even though the smaller recirculation zone appearing [50] near the leeward corner is not predicted here.

Furthermore, instantaneous streamlines that correspond to time instants 1/20, 5/20, 9/20, 11/20, 15/20, and 19/20 of the period (τ) of vortex shedding are shown in Figure 11. These streamlines are in a good agreement with the streamlines computed with the Launder and Sharma low-Re k-ε model by Raisee et al. [66].

3.2. Gust Introduction in Flow Past a Square Cylinder

For implementation of the wind gusts, described in Section 2.3, the inlet velocity profile of Equation (8) and Figure 2 was applied with a period of T₁ = 20 s. and T₂ = 100 s for a wind gust of short and long duration, respectively. These test cases were conducted with free-stream boundary conditions on the two sides in the lateral dimension.

Points a, b, c, d, and e in Figure 2 denote particular instances of time that are denoted in the same manner in Figure 12, Figure 13 and Figure 14. In Figure 12, the contours of the streamwise velocity are presented for the steady-state upstream flow before any gust begins, at position a in Figure 2.

On the other hand, Figure 13b₁,b₂ correspond to time instance b, i.e., at the end of the initial reduction of velocity. Although the absolute decrease in velocity value is the same for the two gusts, the effect of their duration is already evident as, in Figure 13b₁ for the short-duration gust, the absolute value of velocity is smaller in the wake than in Figure 13b₂, where, due to the longer duration, the inlet velocity was reduced at a slower rate.

Figure 13c₁,₂ presents the flow field when the inlet velocity reaches its maximum value. For the longer gust duration (Figure 13c₂), a large region of flow reversal appears in the wake of the rod, which is almost completely absent for the gust of shorter duration (Figure 13c₁). Points b and d in Figure 2 correspond to the same inlet velocity, but the time history of the velocity leads to Figure 13b_1,2,d_1,2 having completely different flow fields.

It is impressive to note the symmetric regions that appear above and below the square cylinder for the shorter-duration gust, which persist until the gust is over (Figure 13e₁). However, when the longer-duration gust is introduced, the rate of change of upstream velocity is slower, and the general form of the natural vortex-shedding sequence is still recognizable. This observation may be examined in relation to the natural vortex-shedding period, which was 23.5 s, slightly longer than the short-duration gust. On the other hand, the gust of longer duration contains more than four natural vortex-shedding periods. This behavior is reminiscent of the lock-on phenomenon, which appears when the frequency of model oscillation, or upstream velocity fluctuation frequency, approaches the natural shedding frequency [67,68].

In Figure 14, contours of streamwise velocity are presented after the end of the gust duration. The letters f, g, h, and i correspond to 25%, 50%, 75%, and 100% of each gust’s period, respectively. Therefore, Figure 14f₁–i₁ and Figure 14f₂–i₂ present results that correspond to less than one natural vortex-shedding period and more than four vortex-shedding periods of each gust, respectively. When the short-duration gust is over, the flow needs more than half of a gust period to recover, and the pair-vortex created on the sides of the rod persists until then.

On the other hand, for the longer-duration gust, the flow recovers almost immediately, and the flow field seems to continue as if the gust never happened. It should be kept in mind though that 25% of the long-duration gust period corresponds to roughly one natural shedding period and/or one short-duration gust period.

In Figure 15, the timeseries of the streamwise velocity component at the centerline, five cylinder side lengths downstream (x/D = 5), is presented. The shape of the gust signal applied at the inlet (Figure 2) is barely recognizable in either of the two curves (the gust begins at t* = 25). While the peak is obvious in both cases, the initial decrease in wind velocity, i.e., at the beginning of the gust, is bigger than the one near the end of the shorter-duration gust, while this is reversed for the longer-duration gust. As was obvious from the contours in Figure 13, velocity fluctuations were larger when the longer-duration gust was introduced. However, the delay of flow recovery is apparent in Figure 15a, after the end of the short-duration gust. Even if one considers that purely periodic behavior cannot be expected in turbulent flows, it is noteworthy that while the basic periodic behavior seems to have recovered at 60 non-dimensional time units after the long-duration gust (Figure 15b), the streamwise velocity is still not periodic at 100-non dimensional time units after the short-duration gust (Figure 15a). It should be noted that the streamwise velocity presented in Figure 15 took longer than all the other variables to recover.

The timeseries of the transverse velocity component at the same above-mentioned position (x/D = 5) is presented in Figure 16. Transverse velocity is significantly reduced, and its recovery is delayed in the case of the short-duration gust, while for the longer duration, the periodicity in the transverse velocity signal immediately recovers.

The total drag coefficient (C_D) was also calculated from the difference between the pressure coefficient on the windward and leeward sides. As the friction coefficient did not have a significant contribution, C_D was calculated from

$${\mathsf{c}}_{\mathsf{D}}=\frac{2}{\mathsf{\rho }\mathsf{A}{\mathsf{u}}_{\mathsf{i}\mathsf{n}\mathsf{l}\mathsf{e}\mathsf{t}}^2}\cdot \left[{\displaystyle \underset{\mathsf{S}}{\int }{\mathsf{P}}_{\mathsf{w}\mathsf{i}\mathsf{n}\mathsf{d}\mathsf{w}\mathsf{a}\mathsf{r}\mathsf{d}}\mathsf{d}\mathsf{s}}-{\displaystyle \underset{\mathsf{S}}{\int }{\mathsf{P}}_{\mathsf{l}\mathsf{e}\mathsf{e}\mathsf{w}\mathsf{a}\mathsf{r}\mathsf{d}}\mathsf{d}\mathsf{s}}\right]$$

(11)

where ρ is air density, A is the cylinder’s frontal area, and u is the mean inlet velocity. The timeseries of total drag coefficient are presented in Figure 17 for the two gusts. For the long-duration gust, the C_D time signal is similar to that of the inlet streamwise velocity.

The maximum C_D value is the same in the two timeseries, but when it comes to the short-duration gust, after reaching the maximum value, its value significantly drops (the corresponding time-instance is that of point d₁ in Figure 2) and eventually reaches a negative value.

The lift coefficient is calculated in a similar way to the drag coefficient through the difference in pressures exerted on the upper and lower side of the cylinder:

$${\mathsf{c}}_{\mathsf{L}}=\frac{2}{\mathsf{\rho }\mathsf{A}{\mathsf{u}}_{\mathsf{i}\mathsf{n}\mathsf{l}\mathsf{e}\mathsf{t}}^2}\cdot \left[{\displaystyle \underset{\mathsf{S}}{\int }{\mathsf{P}}_{\mathsf{u}\mathsf{p}\mathsf{p}\mathsf{e}\mathsf{r}}\mathsf{d}\mathsf{s}}-{\displaystyle \underset{\mathsf{S}}{\int }{\mathsf{P}}_{\mathsf{l}\mathsf{o}\mathsf{w}\mathsf{e}\mathsf{r}}\mathsf{d}\mathsf{s}}\right]$$

(12)

The timeseries of the lift coefficient are presented in Figure 18, and it is notable that throughout the short-duration gust, C_L is almost zero. On the other hand, for the gust of longer duration, the total lift coefficient still maintains its periodicity, but there is a single fluctuation that has a value exceeding those of the previous and following periods. It is evident that the gust influence in the lateral direction is significant and also depends on the duration of the gust, i.e., the frequency of the fluctuation in relation to the natural shedding frequency.

In an effort to further analyze the behavior of the drag and lift coefficients, timeseries of the pressure coefficients on the windward and leeward sides are presented in Figure 19. These correspond to the first and second term of Equation (11), respectively. These curves are expected to be proportional to the C_D curve in Figure 17, but the individual behavior and the phase lag between opposing sides holds important information that is lost in the total drag.

During the short-duration gust, the curves of the upwind and leeward pressure coefficients have a similar shape, but as expected, c_p on the windward side is over twice the value of the leeward side. Interestingly, when the inlet velocity takes its maximum value (according to Figure 2), c_p is significantly increased, as expected, but the subsequent steep drop in inlet velocity corresponds to an extreme drop in c_p, (proportionally bigger than the drop in inlet velocity). On the contrary, when a longer-duration gust occurs, c_p changes smoothly, and although, again, the second steep drop is much bigger than the first one, the values of c_p are approximately half those of the short-duration gust.

The timeseries of pressure coefficients calculated on the upper and lower sides of the square cylinder are presented in Figure 20. Again, these coefficients correspond to the first and second term of Equation (12), respectively, and have mostly negative values due to the suction appearing relative to the free-stream reference pressure. During the short-duration gust, the two coefficients are almost the same, as expected from Figure 18, where the total lift coefficient was zero. Additionally, the values of the two pressure coefficients are significantly reduced when the steep drop in inlet velocity occurs after the gusts peak. However, for a longer-duration gust, C_P changes smoothly on the top and bottom sides, their periodicity is recognizable throughout the gust, and, when it is over, they recover their initial values and periodicity very quickly. Furthermore, the two pressure coefficients continue to have a different phase, contrary to the short-duration gust, where the two pressure coefficients have the same phase and magnitude.

In Figure 21, Lissajous curves are also presented, with C_D on the horizontal axis and C_L on the vertical axis. When the short-duration gust occurs, C_L is zero, and the shape of the Lissajous curve is not maintained during the gust. However, when a longer-duration gust occurs, the velocity slowly changes, and some form of periodicity is maintained even for the duration of the gust. Therefore, the Lissajous curve loop is displaced during the long duration gust, maintaining almost the same shape but returning to its initial position immediately after the gust.

4. Discussion and Conclusions

Single wind gusts of short and long duration were applied upstream of a square cylinder in order to investigate their effect on global characteristics, such as the wake, the drag and lift coefficients, and the pressure distribution on the cylinder’s outer surfaces. Although these characteristics are most often of interest to applications related to the forces acting on the bluff body, when dealing with free-standing constructions such as towers or buildings, the pressure differences are also relevant to infiltration and ventilation.

This research was carried out by implementing a very large eddy simulation (VLES) model to an in-house code, which was then successfully validated against results from the literature for a square cylinder exposed to a steady upstream flow. The validated in-house code was then used to study the effect of two upstream gusts, whose differences lie in their duration but not in the values of the mean inlet velocity or the gust amplitude. The duration of the short gust was very close to the natural vortex-shedding period, and the longer gust duration is four times larger.

The change in behavior of the flow characteristics was found to be significant under the influence of the gusts. Detailed comparison and quantitative discussions have already been presented along with the results, but the applied implications are further discussed here. For the short-duration gust, the behavior closely resembles the lock-on phenomenon frequently visited in case of a periodic upstream velocity fluctuation; i.e., the gust’s sudden change in inlet velocity breaks the alternating vortex-shedding pattern and results in the creation of pair-vortex-shedding, as the two vortices shed simultaneously from the lateral sides. On the other hand, for the longer gust duration, the gradual change in the velocity allows the flow to maintain its vortex-shedding periodicity even though the amplitude of the gust is the same as in the previous case. This conclusion agrees with the expected behavior for upstream perturbation frequencies close (short-duration gust) or far (long-duration gust) to the natural shedding frequency of the cylinder.

It is interesting to discuss the implications of the formation of the pair-vortex during the short gust duration, as it leads to simultaneously occurring, identical pressure coefficients on the lateral sides of the cylinder, and so the overall lift coefficient is close to zero. This occurs even though these pressure coefficients are more than twice as large as those for the gust of long duration. In the case of a hollow structure, such as a building, where air exchange is determined by pressure differences, short-duration gusts would induce minimal cross-ventilation air exchange in the lateral direction.

Conversely, more efficient lateral air exchange is achieved for longer-duration gusts because although the pressure coefficients are of lower magnitude, their variation is out of phase and leads to an overall pressure difference across the building. When it comes to the streamwise direction, the differences between the short- and long-duration gust do not significantly affect the maximum values of the total drag coefficients. However, for the case of the short-duration gust, the total drag coefficient and even the windward pressure coefficient obtains a negative value near the end of the gust, indicating air exfiltration even from the windward side.

An overall observation is that gust occurrence will contribute to air exchange, but depending on the duration of the gust, this may or may not be significant in the streamwise or lateral direction. For the short-duration gust, occurring close to the natural shedding period of the cylinder, streamwise air exchange is favored, notably in both upstream and downstream directions, but lateral air exchange is minimal. Longer-duration gusts may induce significant lateral air exchange and streamwise exchange in the downstream direction. It should be noted that the present study did not examine different gust periods to find the development of the behavior as gust duration increases, nor was the overall exchanged air volume considered, especially for the case of repeated gusts of the same or different duration. These are novel observations that justify further research with relevance to air infiltration and building natural ventilation.