Numerical simulation of vortex vibration in main girder of cable-stayed bridge based on bidirectional fluid–structure coupling

The occurrence of VIV bridges is mainly related to flow separation and vortex shedding, and its fundamental reason lies in the periodic shedding vortex generated when the wind with constant velocity passes through the main girder. When there is a vortex around the structure, it will be subjected to a stable periodic force in the direction perpendicular to the fluid, while when the frequency of the force is close to the inherent frequency of a certain order of the structure, it will cause resonance [11].

When the air flows around the cylinder, it will produce a vortex that periodically falls off at the wake. In this paper, the dimensionless Storrow's number is defined to describe the periodic characteristics of the vortex, namely: (1) St=nDU=f(ρUDμ)=f(Re) St = \frac{nD}{U} = f \left( \frac{ \rho UD}{\mu}\right) = f(\text{Re}) where, n is vortex stripping frequency; D is the characteristic size of the structure perpendicular to the plane of incoming flow direction; U is wind speed of incoming flow; ρ is air density; μ is the coefficient of air viscosity; and Re represents the Reynolds number.

When the air flows through the cross-section of a blunt bridge structure, the wake will have periodic vortex shedding alternately, which will lead to the change in the surface pressure of the bridge and eventually lead to the vibration of the bridge, which is called VIV.

Under the average wind, the bridge structure continuously absorbs energy from the fluid flowing around it to produce vibration, and the vortex that alternately falls off at the wake reacts on the bridge structure to cause vibration, so the VIV is self-excited. While under the wind pulsation, the bridge will cause random vibration due to its pulsating components, it will also have some influence on vortex shedding, thus limiting the amplitude of bridge vortex vibration. Therefore, vortex vibration is a kind of forced vibration with limited amplitude, and it often has more than one wind speed interval. The intensification of vibration leads to strong coupling between cylinder and fluid, and controls VIV in a certain range of wind speed. At this time, the vortex shedding frequency is always consistent with the natural frequency of the cable, as shown in Figure 1.

Fig. 1

Vortex vibration is a common phenomenon of wind-induced vibration in engineering structures, and the quality, damping ratio and natural frequency of the structure will have a greater or lesser impact on its vortex vibration performance. The natural frequency of the structure plays a decisive role in the range of wind speed locked by vortex vibration, while the magnitude of vortex amplitude is affected by the damping ratio and mass of the structure, where the larger the mass of the structure is, the more energy it needs to absorb from the airflow. Meanwhile, the damping ratio determines the energy dissipated when the structure vibrates, so that the larger the damping ratio is, the more energy dissipated [12].

Numerical simulation of bridge vortex vibration based on CFD/CSD can obtain the information of wind load and structural coupling motion comprehensively. Technically, it is necessary to establish CFD flow field and CSD domain model, and consider the interface coupling between flow field and domain. As CFD technology mainly adopts finite difference or finite volume discrete continuous N-S equations, while CSD technology mainly adopts finite element discrete structural dynamic equations, the bidirectional fluid–structure coupling involves the mesh matching and interface movement of the two at the fluid-solid interface.

Bidirectional fluid-solid coupling methods can be divided into strong coupling and weak coupling. In this paper, the vortex vibration of main girder of cable-stayed bridge is numerically simulated by weak coupling, that is, CFD solver and CSD solver solve and exchange data in one time step. The specific process is shown in Figure 2.

First, the surface load result of the current time step model is obtained by CFD solver;

Fig. 2

3.1.2

Simulation steps

In this paper, the Gambit pre-processing software is used to calculate the bidirectional fluid-solid coupling of the steel box girder section of cable-stayed bridge, so as to explore the vortex vibration performance of the girder section. The specific process is as follows:

3.1.2.1

Divide computing domains

In this paper, the Gambit pre-processing software is used to model the steel box girder section of cable-stayed bridge and divide the grid. If the width of the calculation area is not enough, the vortex separated at the tail of the main girder structure model will move to the boundary and then return. Meanwhile, if the length of the calculation area is not enough, it is difficult to accurately define the exit boundary conditions, while too large a calculation area will increase the calculation time.

To simulate the wind environment in the natural state as much as possible, we need to reasonably define the size of the calculation domain and make the calculation basin meet the requirement that the blocking rate is <5%. In this paper, a rectangular computational domain of 28D × 17B is selected, and the computational domain is divided into rigid domain, dynamic grid domain and stable domain. The rigid domain is 1.12D × l.25B in size, while the distances from the main girder to the upper edge, the lower edge, the left edge and the right edge of the rigid body domain are 0.05B, 0.1B, 0.125D and 0.125D, respectively, and the size of the moving grid area is 3B × 10D. To fully develop the wake flow field, the size of the tail basin is selected as 9.5B × 10D.

3.1.2.2

Divide grid

On the premise of ensuring the accuracy of numerical simulation calculation, reducing the number of n grids s as much as possible can shorten a lot of calculation time. In this paper, quadrilateral mesh is used for division. Displacement caused by vibration in rigid area will be transferred to moving grid area, so the grid in moving grid area needs to be reconstructed and updated constantly to ensure that the grid will not be distorted. Therefore, in this paper, the triangular grid is used in the moving area, so that it can be set as the moving grid in the later coupling calculation, and the grid reconstruction method is used to update and reconstruct the grid in the area of moving grid.

Because the wake area is large, to effectively reduce the number of grids, the whole quadrilateral grid is used to divide this area. Then, the boundary layer on the main girder section and ancillary facilities are set to better simulate the fluid situation around the main girder section, so as to obtain more accurate numerical simulation results. At the same time, to verify the grid independence of fluid domain, the whole fluid domain is divided into two sets of grids, which are 182,600 and 360,000, respectively. Through trial calculation, it is found that the grids of 260,000 and 360,000 are in good match with the experimental results, which shows that the increase in the number of grids can no longer affect the calculation results. Therefore, 260,000 grids are finally selected as the final calculation grid where the first grid is 0.016 mm, with 15 floors.

3.1.2.3

Setting of boundary conditions

In this paper, the boundary condition of the entrance is set as the velocity-inlet. In this paper, the lower boundary and the left boundary are set as the entrance boundary, the corresponding right boundary and the upper boundary are set as the exit boundary, and the exit boundary condition is set as the pressure-outlet. In addition, the main girder section and the accessory structure model surface are set as the wall boundary condition without slip. Figure 3 shows the specific boundary conditions.

Fig. 3

Moreover, the steps of boundary setting are shown in Figure 4:

Use SST k-w turbulence model to calculate, and the turbulence intensity and turbulence viscosity ratio are set to 0.5% and 5, respectively;

Adopt SIMPLE algorithm to solve the problem in steady state;

After the steady-state calculation results converge and stabilise, the SIMPLEC algorithm is used for transient calculation;

The velocity-pressure coupling is used to solve the problem, and the second-order upwind scheme is used to control the dispersion of the calculation.

Fig. 4

For the rigid VIV of the steel box girder studied in this paper, the displacement, speed, acceleration, wind speed and time of the box girder will affect the instantaneous pressure acting on the surface of the girder. Without considering the influence of turbulence, VIF can be described by the function of wind speed, time, displacement, speed and acceleration of box girder. In this paper, the rigid VIV of the steel box girder mainly occurs at the main frequency, so the influence of high-order frequency on the displacement can be neglected, where the displacement and acceleration of the box girder have the same fluctuation characteristics except for the difference in phase and amplitude, so a variable can be used to replace these two physical quantities. Therefore, the rigid VIF of the steel box girder can be expressed as follows: (8) fVT=f(y,y˙,U,t) f_{VT} = f(y, \dot{y}, U, t)

The vortex-excitation force consists of four kinds of components: the term containing displacement y only, the term containing velocity ẏ_* only, the cross term of displacement y and velocity ẏ_*, and the vortex-detachment force term. Therefore, the vortex excitation force can be expressed as follows in a broad sense: (9) fVI=ρU2D[∑i=0,j=0i=n,j=nPij(y˙*)i(y˙*)j+Vssin(ωst+ϕ)] f_{VI} = \rho U^{2} D \left[ \sum\limits_{i=0, j=0}^{i=n,j=n} P_{ij} \left( \dot{y}_{*}\right)^{i} \left( \dot{y}_{*}\right)^{j} + V_{s} \sin \left( \omega_{s} t + \phi \right)\right] where ρ represents air density; U represents wind speed of incoming flow; D represents the height of the main beam; y_* = y/D represents the dimensionless vortex displacement, ẏ_* represents the dimensionless vortex vibration velocity (ẏ_* = ẏ/U); P_ij, V_s, ω_s and ϕ represents the undetermined parameters.

It can be seen from Eq. (9) that in addition to the vortex detachment force, the first-order response terms in the model only include the linear aerodynamic stiffness term P₁₀ẏ_* and the linear aerodynamic damping term P₀₁y_*, and the other terms are all nonlinear high-order terms.

To establish the VIF mathematical model of box girder with cantilever more reasonably, it is necessary to accurately judge the influence of terms contained in the generalised expression of VIF. Therefore, in this paper, the residual error between the VIF reconstructed by the mathematical model and the actual VIF is used as the criterion, and its expression is as follows (10) R=1n∑i=0n[f17m(ti)−fV(ti)]2 R = \frac{1}{n} \sum\limits_{i=0}^{n} \left[ f_{17}^{m} (t_{i}) - f_{V}(t_{i}) \right]^{2} where R represents the residual value; n indicates the number of time steps, f17m f_{17}^{m} is the model reconstructed VIF; and f_VI represents the actual VIF. After comprehensive consideration of the items and residuals, the expression of the vortex excitation force mathematical model finally established by the above method is as follows: (11) fV=ρU2D[P10y˙*+P01y*+P11y˙.y*+P12y˙.y*2+P21y˙.y*+Vssin(ωst+ϕ)] f_{V} = \rho U^{2} D \left[ P_{10}\dot{y}_{*} + P_{01} {y}_{*} + P_{11}\dot{y}.y_{*}+ P_{12}\dot{y}.y_{*}^{2} + P_{21}\dot{y}.y_{*}+V_{s}\sin\left( \omega_{s}t + \phi \right)\right]

In this paper, the bidirectional fluid-solid coupling calculation of the two-dimensional main girder section under +3 wind attack angle is carried out by using the same related parameters as the segmental model test, and the numerical simulation results are compared with the segmental model test results. The numerical simulation parameters are shown in Table 2.

The variation of vortex vibration amplitude of the main girder section with wind speed under the action of +3 wind attack angle is as shown in Figure 5.

Fig. 5

Under +3 wind attack angle, the variation trend of vortex amplitude of steel box girder section at different wind speeds obtained by numerical simulation is basically the same as that obtained by the segmental model test, and the vortex vibration locking interval measured by them is basically the same, roughly between 5.5 m/s and 7.5 m/s. On the whole, however, the vortex vibration results obtained by numerical simulation are basically consistent with those measured by two-dimensional segmental model test, whose results are reliable.

The frequency spectrum analysis of numerical simulation results shows that the VIV frequency of numerical simulation is 9.498 Hz, while the natural frequency of the main beam of two-dimensional segment test is 9.745 Hz, and the difference between them is about 2.53%. Therefore, the numerical simulation method adopted in this paper can accurately simulate the VIV process of the steel box main girder, and the simulation results are accurate and reliable.