# Nonlinear Flutter Analysis for a Long-Span Suspension Bridge

Authors: Jieshan Liu; Fan Wang; Yang Yang; Renhuai Liu
Journal: Shock and Vibration
 (2021)

---

## Abstract

The wind-induced flutter of the long-span suspension bridge structure is extremely harmful to the bridge. Therefore, it is necessary to study the nonlinear problem of wind-induced flutter. Here, the nonlinear flutter problem of a long-span suspension bridge with cubic torsional stiffness is analyzed by the equivalent linearization method. The system has Hopf bifurcation and limit cycle oscillations (LCOs) under critical wind speed. Replacing the nonlinear stiffness term of the original nonlinear equation with the equivalent linear stiffness, we can obtain the equivalent linearized equation of the nonlinear flutter system and the solution, critical wind speed, and flutter frequency of the suspension bridge flutter system. At the same time, the system has a limit cycle vibration, and the Hopf bifurcation point is obtained. Compared with the numerical method, the calculation results are consistent. The influence of the damping ratio on the flutter system is analyzed. Increasing the system damping ratio can increase the flutter critical wind speed and reduce the amplitude of LCOs. The influence of cubic torsional stiffness on the flutter system is analyzed. The increase of the cubic stiffness coefficient does not change the critical state of flutter, but reduces the amplitude of LCOs.


---

## Body

## 1. Introduction

With the increase of the bridge span, especially the long-span suspension bridge, the structure tends to soften. Its damping and natural vibration frequency of the structure will decrease. The influence of wind on the bridge can easily cause wind vibration. Among them, flutter is the most serious wind-induced vibration phenomenon because the vibration will radiate when flutter occurs, and the amplitude of vibration will gradually increase with the increase of time and wind speed, which will seriously damage the bridge. Therefore, flutter cannot occur throughout the life cycle of the bridge. The flutter problem can be regarded as a vibration produced by the self-excited forces caused by airflow acting on the bridge deck center, according to Scanlan’s self-excited flutter model; the model assumes that the two-degree bridge deck model is excited by self-excited force when it vibrates [1]. The self-excited force is expressed as a linear function of vibration displacement and velocity in the form of flutter derivative; then, the flutter derivatives are the function of reduced wind speed, and they can be obtained through wind tunnel tests or numerical simulations. The flutter analysis based on linear elasticity and self-excited force theory can accurately calculate the critical state of the bridge, and when the wind speed exceeds the critical flutter wind speed, the bridge vibration will increase rapidly and becomes unstable. This is a structural instability where the equilibrium state of the structure begins to lose its stability under the action of external forces, and a slight disturbance, which is practically unavoidable, gradually increases the deformation and finally causes the structure to break down. The instability of the bridge system at the critical state of flutter is Hopf bifurcation; when wind speed exceeds it, there will be a phenomenon such as the Tacoma Narrows Bridge and the finite-amplitude periodic vibration that occurrs in the wind tunnel test, that is, the limit cycle oscillations (LCOs). LCO is a steady-state motion such as a stable limit cycle that diverges from the static instability point in Hopf bifurcation; in the phase diagram, it is a closed isolated trajectory. Only the critical flutter state and Hopf bifurcation can be calculated if based on a linear flutter system. However, the LCO is the solution of a nonlinear system. Therefore, a nonlinear system based on the nonlinear characteristics of large-span suspension bridges is needed to analyze the Hopf bifurcation type and LCO after the bifurcation.On the one hand, for many long-span suspension bridges, the effects of geometric and structural nonlinearities are prominent due to the reduced structural stiffness, making the bridge section vibrations negligible. On the other hand, structurally, the steel box girder in the bridge system is similar to the aircraft wing as in Figure1, in which Figure 1(a) is a construction of the steel box girder and Figure 1(b) is a construction of the aircraft wing. The nonlinear flutter problem of long-span suspension bridges is analogous to the flutter problems of large aspect ratio flexible wings, which should be studied theoretically using the nonlinear dynamics analysis method for the Hopf bifurcation problem, in which the limit cycle oscillations (LCOs) occur after the bridge system enters into a flutter.Figure 1
(a) Structure of the steel box girder of the bridge system. (b) Structure of the aircraft wing.
(a)(b)At present, the nonlinear factors of nonlinear flutter can be divided into three categories: aerodynamic self-excited force nonlinearity, geometric nonlinearity, and structural nonlinearity. Flutter was initially studied in the field of aircraft. After the Tacoma Narrow Bridge broke in the wind in 1940, experts used aircraft theory to study the flutter of the bridge and then developed it into a bridge aerodynamic theory.In the study of nonlinear aerodynamic self-excited forces, Scanlan introduced the self-excited forces, which can be expressed as a function of the vertical and torsion velocity and displacement of the bridge deck, and the coefficients are expressed as the flutter derivatives which associated with reduced wind speed [1–5]. Scanlan pointed out that the aerodynamic self-excited forces of the bridge deck are not only a function of reduced wind speed but also related to the wind attack angle and turbulence of the bridge deck [6]. The results of the wind tunnel test and numerical simulation show that aerodynamic self-excited force has a great relationship with wind attack angle and turbulence [7–12]. Also, the nonlinearity of self-excited aero force can be analyzed by the numerical finite element method [13, 14].In the study of nonlinear problems caused by the geometric factor of suspension bridges, it is a system composed of two-degree bridge deck and spring constant hangers. Mckenna established this geometrically nonlinear system [15]. When the bridge deck is subjected to torsional displacement, the change in the magnitude and direction of the elastic force of the hanger leads to geometric nonlinearity of torsion. Based on this system, the authors studied the bifurcation and LCOs of this kind of two-degree bridge deck model [15–21].When analyzing the nonlinear problems caused by the nonlinear stiffness of the structure, the first step is to study the aircraft. Woolston once said that when the wing twists, its stiffness coefficient is a nonlinear function of the torsional angle, which is called cubic, flat spot, or hysteresis nonlinear stiffness springs [22]. In addition, Woolston also made a spring mechanism for the wings, which can be adjusted to show three types of nonlinear stiffness characteristics in order to measure its free vibration performance in a wind tunnel, and then discussed the relationship between wind speed and amplitude [23]. The cubic nonlinear stiffness coefficient is written as a cubic function of the torsion angle and has been widely used in the flutter equations of aircraft [24–29].Dowell concluded that the nonlinear stiffness wing will exhibit LCOs when wind speed exceeds the critical wind speed in many wind tunnel tests and numerical simulations [30]. Wright recapitulated that structural nonlinearities will occur when the aircraft’s components undergo large deformations, and they will exhibit cubic nonlinear stiffness when the wings undergo large deformations [31]. Chen discussed the problem of flutter LCOs of a wing with cubic nonlinear stiffness by ELM [27]. Chen established a torsional stiffness term with cubic nonlinear stiffness by determining the characteristic value of the linear system to obtain the critical flutter wind speed, and the LCO of a wind speed value was solved by the iteration method [32].The bridge deck system of the long-span suspension bridge mostly adopts the steel box girder structure. To some extent, the steel box girder structure is similar to the large aspect ratio wing structure of aircraft. Both of them have longitudinal beams with web plate and diaphragms as in Figure1. The webs contribute greatly to the structural rigidity. According to Younger’s web shear force distribution theory, the shear stress of the web presents a nonlinear distribution when H-shaped and I-shaped cross-section beams deform in bending and torsion [33, 34]. On this basis, Kuhn Paul studied the longitudinal rod-web combined model. The longitudinal force acts on the longitudinal rod, so the longitudinal rod transmits force to the web, and the reaction force of the web affected it. Therefore, the relationship between force and displacement is obtained in line with Younger’s theory. The longitudinal force and torsion force will show nonlinear changes with the increase of displacement [35]. Ma conducted a stiffness test on the deck system of the suspension bridge, and its torsional stiffness increased with the torsion angle during the small displacement stage, showing a strong nonlinear growth trend [36]. Therefore, this article holds that the webs in the steel box girder have a nonlinear stiffness enhancement effect on the overall structure so that the steel box girder exhibits nonlinear torsional stiffness.In this paper, the flutter and LCOs beyond critical flutter speed of a long-span suspension bridge deck with cubic nonlinear torsional stiffness will be studied. First, based on Scanlan’s flutter self-excited forces theory, the discrete test flutter derivatives were fitted to the quadratic function of the reduced wind speed. Second, considering the cubic torsional stiffness and the vibration frequency of a long-span suspension bridge, nonlinear flutter equations of Nizhou Waterway Bridge were established by referring the plane wings flutter model. Third, by replacing the nonlinear stiffness term with the equivalent linear coefficient through ELM, the critical flutter wind speed and critical flutter frequency are obtained, and at the same time, at the point of critical flutter wind speed, Hopf bifurcation and LCOs occur.


## 2. Nonlinear Flutter Model of Bridge Deck

The bridge flutter model is assumed as a two-dimensional rigid bridge deck suspended by two springs and oscillated in vertical (denoted byh) and torsional (denoted by α) displacements, as shown in Figure 2. Two springs have vertical stiffness kh and torsional stiffness kα, respectively; thus, the vertical bending stiffness and torsional stiffness of a large-span suspension bridge are expressed by each springs’ stiffness. The stiffness coefficients of the springs have a relationship with the first-order frequency ωh or ωα and mass mh or I of the actual bridge in vertical bending and torsion, that is, kh=mhωh2 or kα=Iωα2. Each bridge has different vibration frequencies, especially for a large-span suspension bridge, where the first-order frequencies in vertical and torsional are smaller than for other bridge types. Therefore, in this paper, the span of 1688 m, named Nizhou Waterway Bridge, was selected as the object of study; the stiffness coefficients in the flutter model are related to the first-order self-oscillation frequency of the whole bridge. According to Scanlan’s flutter theory [1] and research in aircraft [37, 38], the self-excited forces of two-degree (vertical excited force denoted by FL, torsional excited force denoted by FM) are acting on bridge deck center.Figure 2
Two-degree bridge deck model.This bridge flutter model is widely used in many bridge flutter studies [4, 5, 9, 39–41]; it was altered from the aircraft wing flutter model by many scholars in the past.


### 2.1. The Cubic Nonlinear Stiffness

The bridge deck will show structural nonlinearity when taking large torsional amplitude, which can be divided into centralized and distributed nonlinearity. In this paper, the centralized nonlinearity in the structural nonlinearity of the long-span suspension bridge has been studied.When analyzing the concentrated nonlinear stiffness of the structure, Ko Jeonghwan and Strganac Thomas expressed the nonlinear force as the force-deformation curve fitting and took the appropriate term. This article uses a concentrated nonlinear form for calculation and research [42, 43].When the bridge deck has a torsional deformation, the torsional stiffness term will be related to the cubic of torsional displacement, which can also describe that it is the cubic nonlinearity stiffnesskα=Iωα2+eαα3. It was represented as the torsional spring having cubic stiffness acting on the flutter model in Figure 2.The flutter motion system of the bridge deck with cubic nonlinear stiffness is written as equations (1) and (2). On the left side of the equations are two-dimensional effect. On the right side of the equations are the flutter self-excited forces which are put forward by Scanlan [1]. The torsional stiffness term on the left side of the equations is written as a cubic function of torsional displacement, which means that the nonlinear cubic stiffness is related to torsional displacement:(1)mhh¨+2ξhmhωhh˙+mhωh2h=ρU2BKH1∗h˙U+KH2∗Bα˙U+K2H3∗α+K2H4∗hB,(2)Iα¨+2ξαIωαα˙+Iωα2α+eαα3=ρU2B2KA1∗h˙U+KA2∗Bα˙U+K2A3∗α+K2A4∗hB,where h, h˙, and h¨ are the displacement, velocity, and acceleration of the vertical direction, respectively, α, α˙, and α¨ are the displacement, velocity, and acceleration of the torsional direction, ξh and ξα are integral damping ratios, ωh and ωα are the corresponding natural mechanical frequencies of vertical and torsional, Hi∗ and Ai∗i=1,2,3,and4 are coefficients that provided for the self-excited aerodynamic force in a linear manner, K is the reduced frequency, and K=Bω/U; in addition, the reduced wind speed is Vr=2π/K=U/Bf, f is vibration frequency, ρ, U, and B stand for air density, wind speed, and width of the bridge deck, and eα is the coefficient of cubic stiffness.In this paper, Nizhou Waterway Bridge of Second Humen Bridge, a long-span suspension bridge with 1688 meters main span, has been taken as the research object. The correlation parameters [32, 44] are shown in Table 1.Table 1
Parameters of the bridge flutter system.
NameUnitValuemkg/m65649Ikg×m7808150Bm41.7ξh(None)0.005ξα(None)0.005ωhHz0.07241ωαHz0.21213eα(None)1000mhωα2


### 2.2. Fitting Flutter Derivatives

When a suspension bridge is suffered from wind flow, the section of the main girder is subjected to aerodynamic forces in the flow field, which are aerodynamic liftL and aerodynamic torque M, and they are expressed in the plural form of aerodynamic forces:(3)L=ρω2b3RICLeiωt,M=ρω2b4RICMeiωt.In equation (3), RI denotes taking the real part in square brackets, where CL and CM are the peak values of aerodynamic force, which are a function of the reduced frequency or reduced wind speed for nonconstant flows, and they can be expressed by the harmonic aerodynamic coefficients Lα, Lh, Mα, and Mh:(4)CL=πLαα+Lhhb,CM=πMαα+Mhhb.In the Scanlan aerodynamic model of self-excited force, flutter derivatives are used to represent the harmonic aerodynamic coefficients, so the flutter derivative should also be a function of the reduced wind speed (denoted byVr) or reduced frequency (denoted by K); they are only associated to the shape of the cross section of the bridge girder. The values of flutter derivatives for reduced wind speed are identified by the wind tunnel or numerical simulation, and they are discrete points; thus, to solve the flutter system in this paper, these discrete points are fitted to the flutter derivative as a function of reduced wind speed by the least squares method:(5)Ai∗,Hi∗=fVr=fK,where Ai∗ and Hi∗i=1,2,3,and4 represent eight flutter derivatives. In this paper, the discrete flutter derivatives, from Hong’s numerical test with the force vibration method [44], were fitted for polynomial functions of reduced wind speed:(6)Hi∗=fVr=Hi1∗×Vr+Hi2∗×Vr2,(7)Ai∗=fVr=Ai1∗×Vr+Ai2∗×Vr2,where Hi1∗, Hi2∗, Ai1∗, and Ai2∗ are coefficients of functions (6) and (7), as shown in Table 2.Table 2
Coefficients of flutter derivatives’ function.
Flutter derivativesHi1∗ or Ai1∗Hi2∗ or Ai2∗H1∗−0.2583−0.005808H2∗−0.20320.01521H3∗0.07225−0.06389H4∗0.2041−0.02397A1∗0.056530.00007511A2∗−0.013230.0002376A3∗0.023690.006157A4∗−0.010450.002351The fitting curves of eight flutter derivatives are shown in Figure3, in which the dots are discrete values via the test and solid lines are fitting curves.Figure 3
Fitting curves of flutter derivatives.
(a)(b)(c)(d)(e)(f)(g)(h)


## 2.1. The Cubic Nonlinear Stiffness

The bridge deck will show structural nonlinearity when taking large torsional amplitude, which can be divided into centralized and distributed nonlinearity. In this paper, the centralized nonlinearity in the structural nonlinearity of the long-span suspension bridge has been studied.When analyzing the concentrated nonlinear stiffness of the structure, Ko Jeonghwan and Strganac Thomas expressed the nonlinear force as the force-deformation curve fitting and took the appropriate term. This article uses a concentrated nonlinear form for calculation and research [42, 43].When the bridge deck has a torsional deformation, the torsional stiffness term will be related to the cubic of torsional displacement, which can also describe that it is the cubic nonlinearity stiffnesskα=Iωα2+eαα3. It was represented as the torsional spring having cubic stiffness acting on the flutter model in Figure 2.The flutter motion system of the bridge deck with cubic nonlinear stiffness is written as equations (1) and (2). On the left side of the equations are two-dimensional effect. On the right side of the equations are the flutter self-excited forces which are put forward by Scanlan [1]. The torsional stiffness term on the left side of the equations is written as a cubic function of torsional displacement, which means that the nonlinear cubic stiffness is related to torsional displacement:(1)mhh¨+2ξhmhωhh˙+mhωh2h=ρU2BKH1∗h˙U+KH2∗Bα˙U+K2H3∗α+K2H4∗hB,(2)Iα¨+2ξαIωαα˙+Iωα2α+eαα3=ρU2B2KA1∗h˙U+KA2∗Bα˙U+K2A3∗α+K2A4∗hB,where h, h˙, and h¨ are the displacement, velocity, and acceleration of the vertical direction, respectively, α, α˙, and α¨ are the displacement, velocity, and acceleration of the torsional direction, ξh and ξα are integral damping ratios, ωh and ωα are the corresponding natural mechanical frequencies of vertical and torsional, Hi∗ and Ai∗i=1,2,3,and4 are coefficients that provided for the self-excited aerodynamic force in a linear manner, K is the reduced frequency, and K=Bω/U; in addition, the reduced wind speed is Vr=2π/K=U/Bf, f is vibration frequency, ρ, U, and B stand for air density, wind speed, and width of the bridge deck, and eα is the coefficient of cubic stiffness.In this paper, Nizhou Waterway Bridge of Second Humen Bridge, a long-span suspension bridge with 1688 meters main span, has been taken as the research object. The correlation parameters [32, 44] are shown in Table 1.Table 1
Parameters of the bridge flutter system.
NameUnitValuemkg/m65649Ikg×m7808150Bm41.7ξh(None)0.005ξα(None)0.005ωhHz0.07241ωαHz0.21213eα(None)1000mhωα2


## 2.2. Fitting Flutter Derivatives

When a suspension bridge is suffered from wind flow, the section of the main girder is subjected to aerodynamic forces in the flow field, which are aerodynamic liftL and aerodynamic torque M, and they are expressed in the plural form of aerodynamic forces:(3)L=ρω2b3RICLeiωt,M=ρω2b4RICMeiωt.In equation (3), RI denotes taking the real part in square brackets, where CL and CM are the peak values of aerodynamic force, which are a function of the reduced frequency or reduced wind speed for nonconstant flows, and they can be expressed by the harmonic aerodynamic coefficients Lα, Lh, Mα, and Mh:(4)CL=πLαα+Lhhb,CM=πMαα+Mhhb.In the Scanlan aerodynamic model of self-excited force, flutter derivatives are used to represent the harmonic aerodynamic coefficients, so the flutter derivative should also be a function of the reduced wind speed (denoted byVr) or reduced frequency (denoted by K); they are only associated to the shape of the cross section of the bridge girder. The values of flutter derivatives for reduced wind speed are identified by the wind tunnel or numerical simulation, and they are discrete points; thus, to solve the flutter system in this paper, these discrete points are fitted to the flutter derivative as a function of reduced wind speed by the least squares method:(5)Ai∗,Hi∗=fVr=fK,where Ai∗ and Hi∗i=1,2,3,and4 represent eight flutter derivatives. In this paper, the discrete flutter derivatives, from Hong’s numerical test with the force vibration method [44], were fitted for polynomial functions of reduced wind speed:(6)Hi∗=fVr=Hi1∗×Vr+Hi2∗×Vr2,(7)Ai∗=fVr=Ai1∗×Vr+Ai2∗×Vr2,where Hi1∗, Hi2∗, Ai1∗, and Ai2∗ are coefficients of functions (6) and (7), as shown in Table 2.Table 2
Coefficients of flutter derivatives’ function.
Flutter derivativesHi1∗ or Ai1∗Hi2∗ or Ai2∗H1∗−0.2583−0.005808H2∗−0.20320.01521H3∗0.07225−0.06389H4∗0.2041−0.02397A1∗0.056530.00007511A2∗−0.013230.0002376A3∗0.023690.006157A4∗−0.010450.002351The fitting curves of eight flutter derivatives are shown in Figure3, in which the dots are discrete values via the test and solid lines are fitting curves.Figure 3
Fitting curves of flutter derivatives.
(a)(b)(c)(d)(e)(f)(g)(h)


## 3. The ELM and the Critical Flutter State

### 3.1. The Equivalent Linearization Method (ELM)

The equivalent linearization method (ELM) is a nonlinear solution based on the idea of the averaging method. The nonlinear system shows the characteristics of approximate harmonics so, here, it mainly analyzes the first-order harmonics in the system motion and transforms the nonlinear system into an equivalent linear system composed of various equivalent parameters (mass, damping, and stiffness). Then, the linear solution method was used to discuss the nonlinear system.According to the calculation method [45], nonlinear equation (2) of the flutter system is written as follows:(8)Iα¨+Iωα2α+fα˙,α=0,where fα,α˙ is the sum of all terms except inertial force and linear elastic force; it includes the damping and cubic nonlinear stiffness in the nonlinear equation, and here, it contains nonlinear damping, but it is not considered in this article. Then, expand fα,α˙ into Fourier series:(9)fα˙,α=r0+∑n=1∞sncosnφ+tnsinnφ,n=1,2,3,….When fluttering, the ELM assumes that the first-order harmonic force has made a majority in vibration, therefore the first-order harmonic force is retained. The first-order coefficients of Fourier series expansion are as follows:(10)r0=12π∫02πfα˙,αdφ,s1=1π∫02πfα˙,αcosφdφ,t1=1π∫02πfα˙,αsinφdφ.The nonlinear equation can be denoted with the equivalent linearized equation as follows:(11)Ieα¨+ceα˙+keα=0.And, suppose that the solution of equivalent linearized equation (11) is(12)α=Acosωet+θ=Acosφ.Among them,ωe is the vibration frequency and A is the amplitude of equivalent linearized equation (11). Hypothetical solution (12) and first-order Fourier series (9) of fα˙,α are substituted into nonlinear equation (8); then, equation (8) is organized into the form of equivalent linear equation (11) and written in the following equation:(13)Iα¨+−1πaωe∫02πfα˙,αsinφdφα˙+Iωα2+1πa∫02πfa˙,αcosφdφα=0.In this article, damping is expressed linearly so that the equivalent calculation of cubic nonlinear stiffness is adopted. And, the equivalent linear stiffness coefficient is obtained by the following calculation method:(14)kαeq=ke=Iωα2+1πa∫02πfα,α˙cosφdφ.Among them,fα˙,α contains the cubic nonlinear stiffness. Assuming that when flutter occurs, the simple harmonic vibration amplitude of the torsional freedom is A, so the equivalent stiffness coefficient is calculated as follows:(15)kαeq=Iωα2+34eαA2,where kαeq becomes the function of A which is the amplitude of the torsional degree.For the convenience of calculation and expression, let(16)kheq=mhωh2.The nonlinear torsional stiffnessIωα2α+eαα3 and vertical stiffness mhωh2h in flutter equations (1) and (2) are substituted by the equivalent linearized stiffness coefficients kαeqα and kheqh:(17)mhh¨+2ξhmhωhh˙+kheqh=ρU2BKH1∗h˙U+KH2∗Bα˙U+K2H3∗α+K2H4∗hB,(18)Iα¨+2ξαIωαα˙+kαeqα=ρU2B2KA1∗h˙U+KA2∗Bα˙U+K2A3∗α+K2A4∗hB.These are the equivalent linearized equations.


### 3.2. The Critical Flutter State

Assuming that when the nonlinear flutter system bifurcates in the critical state, the system vibrates in harmonic motion, and its limit cycle oscillation solution can be written as follows:(19)h=Heiωt,α=Aeiωt,where H and A are amplitudes of the vertical and torsional degree, respectively.Then, these defined solutions are substituted into equivalent linearized equations (17) and (18), then expand reduced frequency K=Bω/U, and divide eiωt:(20)−ω2Hmh+2ξhmhωhωHi+kheqH=ρB2H1∗Hω2i+ρB3H2∗Aω2i+ρB3H3∗Aω2+ρB2H4∗Hω2,−ω2AI+2ξαIωαωAi+kαeqA=ρB3A1∗Hω2i+ρB4A2∗Aω2i+ρB4A3∗Aω2+ρB3A4∗Hω2.In addition, packing like terms,H and A, these equations can be written as homogeneous equations:(21)−ω2mh+2ξhmhωhωi+kheq−ρB2H1∗ω2i−ρB2H4∗ω2H+−ρB3H2∗ω2i−ρB3H3∗ω2A=0,(22)−ρB3A1∗ω2i−ρB3A4∗ω2H+−ω2I+2ξαIωαωi+kαeq−ρB4A2∗ω2i−ρB4A3∗ω2A=0,which can be written as matrix multiplication simplified:(23)H1A1H2A2HA=0,in which H1A1H2A2 is(24)−ω2mh+2ξhmhωhωi+kheq−ρB2H1∗ω2i−ρB2H4∗ω2−ρB3H2∗ω2i−ρB3H3∗ω2−ρB3A1∗ω2i−ρB3A4∗ω2−ω2I+2ξαIωαωi+kαeq−ρB4A2∗ω2i−ρB4A3∗ω2.When the coefficient determinantH1A1H2A2=0, there is a nonzero solution to the homogeneous equation system, and the amplitudes H and A are not zero at this time. It is considered that the bridge system begins to flutter.So, the coefficient determinant is expanded into an equation forω that includes real and imaginary numbers:(25)Fω=fRω+fIωi=0.The real and imaginary terms in expansion equation (25) are arranged into the real part equation fRω and imaginary part equation fIω:(26)fRω=R4ω4+R3ω3+R2ω2+R1ω+R0=0,(27)fIω=I3ω3+I2ω2+I1ω+I0=0.These two equations are quartic equations, and their roots have to beω≠0. Therefore, the imaginary part equation has no constant term and can be written as cubic equation (27).The coefficients of real and imaginary part equations are written as follows:(28a)R4=Imh+ρIB2H4∗+ρmhB4A3∗+ρ2B6A3∗H4∗−ρ2B6A4∗H3∗+ρ2B6A1∗H2∗−ρ2B6A2∗H1∗,(28b)R3=2ρB2IξαωαH1∗+2ρB4mhξhωhA2∗,(28c)R2=−Ikheq−mhkαeq−ρB4kheqA3∗−ρB2kαeqH4∗−4Imhξhξαωhωα,(28d)R1=0,(28e)R0=kheqkαeq,(28f)I3=ρIB2H1∗+ρmhB4A2∗−ρ2B6A1∗H3∗+ρ2B6A3∗H1∗+ρ2B6A2∗H4∗−ρ2B6A4∗H2∗,(28g)I2=−2mhIωhξh−2mhIωαξα−2ρIB4ωαξαH4∗−2ρmhB4ωhξhA3∗,(28h)I1=ρB4kheqA2∗−ρB2kαeqH1∗,(28i)I0=2kheqIωαξα+2kαeqmhωhξh.In this section, the method to solve equivalent linearized equations (17) and (18) is discussed in which the equivalent stiffness coefficients kheq and kαeq are assumed as constant temporarily.In coefficient matrix (24), it contains the unknown ω and eight flutter derivatives Hi∗ and Ai∗i=1,2,3,and4. Since flutter derivatives Hi∗ and Ai∗ are functions of reduced wind speed Vr, expanding equation (25) or real and imaginary part equations (26) and (27) can be regarded as functions of the independent variable ω which includes parameter Vr.In order to find the solution of equations (26) and (27), the following calculation process is required.To do a cycle solving progress, first, a reduced wind speed value is set to obtain eight flutter derivatives and make the real and imaginary part equations onlyω unknown. Then, the positive roots ωR and ωI can be obtained by solving the two equations, respectively. By increasing the reduced wind speed value from zero gradually, several sets of solutions ωR and ωI related to the reduced wind speed can be obtained. Next, the solutions are drawn on the grid in Figure 4; the value of reduced wind speed as abscissa and the value of ωR and ωI as ordinate. Then, all the dots of ωR and ωI are connected to be two lines, and they intersected at a point, which means that the real and imaginary part equations have a same solution; the frequency is ωC=ωR=ωI, and the reduced wind speed is VCr=Vr. So, when Vr=VCr, ω=ωC is the solution of expanding equation (25). Finally, Vcr and ωc are substituted into the coefficient matrix of homogeneous equations (21) and (22); the ratio of vertical to torsional amplitude, H/A, is obtained. In this case, the coefficient determinant equals to zero, and the solutions H and A of homogeneous equations are not equal to zero. It was the solution of the equivalent linearized flutter equations.Figure 4
The intersection of the real part root and imaginary part root.About the intersection point, the value of abscissa is the critical flutter reduced wind speedVCr=10.088, and the value of ordinate is the critical flutter vibration frequency ωC=0.178. Then, the critical flutter wind speed has been calculated that U=75.114 according to Vr=U/Bω, and the vertical-torsional ratio is H/A=9.443, that is, the critical state of the flutter system in this paper. In other words, when the wind speed is 75.114m/s, the suspension bridge will flutter and be instable to bifurcate, and its frequency of vibration is ω=0.178, and its torsional-vertical ratio is H/A=9.443.


### 3.3. Solution of Equivalent Linearization Equations

In this section, the equivalent stiffness coefficientskheq and kαeq in equations (15) and (16) will be taken into account in the solution process because the coefficient kαeq is the function of torsional amplitude A.First, a torsional amplitudeA is the limit cycle amplitude of the system when the flutter occurs is assumed, so the equivalent stiffness coefficient kαeq is determined; the following calculation steps are the same as in Section3.2.The critical flutter solutions including reduced wind speedVCr, frequency ωC, wind speed UC, vertical-torsional ration H/A, and vertical amplitude H can be all obtained. So, these solutions can be interpreted as that the nonlinear flutter system will cause limit cycle oscillation (LCO) when the wind speed is UC. The vertical and torsional amplitude of limit cycle oscillation are H and A, respectively, and the frequency of vibration is ωC.Based on the above analysis and calculation, the characteristic curves of the nonlinear flutter system are drawn in Figures5–8. In order to validate the equivalent linearization method (ELM), numerical solutions by using the Runge–Kutta method [46, 47] are obtained, whose results are shown in Figures 5–8 too.Figure 5
Curve of wind speed and amplitudes of LCOs.Figure 6
Curve of wind speed and frequency of LCOs.Figure 7
Curve of frequency and amplitudes of LCOs.Figure 8
Curve of reduced wind speed and vertical-torsional ratio.As shown in Figure5, the values of wind speed U are taken as abscissa, and the values of vertical and torsional amplitudes are taken as ordinate; then, the same, as in the abovementioned paragraph, two curves in the grid can be obtained. The nonlinear system has flutter instability when critical flutter wind speed is UC=75.114m/s and system bifurcates to LCOs at the same time. Beyond the critical flutter wind speed, the amplitude of LCOs in vertical and torsional degrees increases monotonously when the wind speed is increasing. Also, the solution by ELM agrees well with the numerical results by the Runge–Kutta method.As shown in Figure6, the values of wind speed U are taken as abscissa, and the values of flutter LCOs’ frequency ω are taken as ordinate, and the nonlinear system reaches the critical state of flutter at wind speed UC=75.114m/s, and the critical frequency of flutter is ωC=0.1786. Then, the limit cycle frequency increases monotonically with the wind speed after the bifurcation. Compared with the numerical solution of the Runge–Kutta method at the “×” point in the figure, the frequency of the limit cycle oscillation obtained by the ELM is the larger one.As shown in Figure7, the values of flutter LCOs’ frequency ω are taken as abscissa, and the values of vertical and torsional amplitudes are taken as ordinate; then, two curves are lined through dots of the vertical and torsional degree, respectively, in the grid. The frequency of the LCO of the nonlinear system at the bifurcation point is the critical frequency ωC=0.178. Beyond the critical flutter value, the amplitude of LCOs in vertical and torsional degrees increases monotonously with increasing frequency of LCOs. And, it is that the numerical results are slightly larger than the solutions by ELM.As shown in Figure8, the values of reduced wind speed are taken as abscissa, and the values of the vertical-torsional ratio are taken as ordinate; one curve is lined through the solid square dots. The nonlinear system has flutter instability when critical flutter reduced wind speed VCr=10.088. The vertical-torsional ratio of LCOs decreases monotonously with reduced wind speed.Finally, to validate the feasibility of ELM, the vertical amplitude, torsional amplitudes, and frequency, which at the reduced wind speedVr=10.248, are substituted into the assumed solution h=Heiωt and α=Aeiωt in Section 3.3. Now, the LCO can be drawn into the phaser diagrams shown in Figure 9. The abscissa is the displacement of vertical and torsional degrees, and the ordinate is velocity of vertical and torsional degrees. The solid lines are solutions of ELM, and dots are numerical results. The comparisons show that the results of ELM and numerical method are very closed in the displacement phase, which indicated that the first-order harmonic is dominant in LCO and that the ELM considering the first-order harmonic is practicable. The LCO, which is assumed to be harmonic vibration, is very close to the numerical results. Compared with the LCO of the vertical and torsional degree, the LCO of the vertical degree is closer to the numerical results in displacement and velocity, and the velocity in LCO of the torsional degree is larger than the numerical results with an error of 16% which may be related to the removal of higher-order harmonic vibration in ELM.Figure 9
Limit cycle oscillation (LCO) vs. numerical solutions. (a) LCO in the vertical degree. (b) LCO in the torsional degree.
(a)(b)Therefore, the flutter bifurcation problem of the nonlinear system can be analyzed by the calculation result of the ELM.


## 3.1. The Equivalent Linearization Method (ELM)

The equivalent linearization method (ELM) is a nonlinear solution based on the idea of the averaging method. The nonlinear system shows the characteristics of approximate harmonics so, here, it mainly analyzes the first-order harmonics in the system motion and transforms the nonlinear system into an equivalent linear system composed of various equivalent parameters (mass, damping, and stiffness). Then, the linear solution method was used to discuss the nonlinear system.According to the calculation method [45], nonlinear equation (2) of the flutter system is written as follows:(8)Iα¨+Iωα2α+fα˙,α=0,where fα,α˙ is the sum of all terms except inertial force and linear elastic force; it includes the damping and cubic nonlinear stiffness in the nonlinear equation, and here, it contains nonlinear damping, but it is not considered in this article. Then, expand fα,α˙ into Fourier series:(9)fα˙,α=r0+∑n=1∞sncosnφ+tnsinnφ,n=1,2,3,….When fluttering, the ELM assumes that the first-order harmonic force has made a majority in vibration, therefore the first-order harmonic force is retained. The first-order coefficients of Fourier series expansion are as follows:(10)r0=12π∫02πfα˙,αdφ,s1=1π∫02πfα˙,αcosφdφ,t1=1π∫02πfα˙,αsinφdφ.The nonlinear equation can be denoted with the equivalent linearized equation as follows:(11)Ieα¨+ceα˙+keα=0.And, suppose that the solution of equivalent linearized equation (11) is(12)α=Acosωet+θ=Acosφ.Among them,ωe is the vibration frequency and A is the amplitude of equivalent linearized equation (11). Hypothetical solution (12) and first-order Fourier series (9) of fα˙,α are substituted into nonlinear equation (8); then, equation (8) is organized into the form of equivalent linear equation (11) and written in the following equation:(13)Iα¨+−1πaωe∫02πfα˙,αsinφdφα˙+Iωα2+1πa∫02πfa˙,αcosφdφα=0.In this article, damping is expressed linearly so that the equivalent calculation of cubic nonlinear stiffness is adopted. And, the equivalent linear stiffness coefficient is obtained by the following calculation method:(14)kαeq=ke=Iωα2+1πa∫02πfα,α˙cosφdφ.Among them,fα˙,α contains the cubic nonlinear stiffness. Assuming that when flutter occurs, the simple harmonic vibration amplitude of the torsional freedom is A, so the equivalent stiffness coefficient is calculated as follows:(15)kαeq=Iωα2+34eαA2,where kαeq becomes the function of A which is the amplitude of the torsional degree.For the convenience of calculation and expression, let(16)kheq=mhωh2.The nonlinear torsional stiffnessIωα2α+eαα3 and vertical stiffness mhωh2h in flutter equations (1) and (2) are substituted by the equivalent linearized stiffness coefficients kαeqα and kheqh:(17)mhh¨+2ξhmhωhh˙+kheqh=ρU2BKH1∗h˙U+KH2∗Bα˙U+K2H3∗α+K2H4∗hB,(18)Iα¨+2ξαIωαα˙+kαeqα=ρU2B2KA1∗h˙U+KA2∗Bα˙U+K2A3∗α+K2A4∗hB.These are the equivalent linearized equations.


## 3.2. The Critical Flutter State

Assuming that when the nonlinear flutter system bifurcates in the critical state, the system vibrates in harmonic motion, and its limit cycle oscillation solution can be written as follows:(19)h=Heiωt,α=Aeiωt,where H and A are amplitudes of the vertical and torsional degree, respectively.Then, these defined solutions are substituted into equivalent linearized equations (17) and (18), then expand reduced frequency K=Bω/U, and divide eiωt:(20)−ω2Hmh+2ξhmhωhωHi+kheqH=ρB2H1∗Hω2i+ρB3H2∗Aω2i+ρB3H3∗Aω2+ρB2H4∗Hω2,−ω2AI+2ξαIωαωAi+kαeqA=ρB3A1∗Hω2i+ρB4A2∗Aω2i+ρB4A3∗Aω2+ρB3A4∗Hω2.In addition, packing like terms,H and A, these equations can be written as homogeneous equations:(21)−ω2mh+2ξhmhωhωi+kheq−ρB2H1∗ω2i−ρB2H4∗ω2H+−ρB3H2∗ω2i−ρB3H3∗ω2A=0,(22)−ρB3A1∗ω2i−ρB3A4∗ω2H+−ω2I+2ξαIωαωi+kαeq−ρB4A2∗ω2i−ρB4A3∗ω2A=0,which can be written as matrix multiplication simplified:(23)H1A1H2A2HA=0,in which H1A1H2A2 is(24)−ω2mh+2ξhmhωhωi+kheq−ρB2H1∗ω2i−ρB2H4∗ω2−ρB3H2∗ω2i−ρB3H3∗ω2−ρB3A1∗ω2i−ρB3A4∗ω2−ω2I+2ξαIωαωi+kαeq−ρB4A2∗ω2i−ρB4A3∗ω2.When the coefficient determinantH1A1H2A2=0, there is a nonzero solution to the homogeneous equation system, and the amplitudes H and A are not zero at this time. It is considered that the bridge system begins to flutter.So, the coefficient determinant is expanded into an equation forω that includes real and imaginary numbers:(25)Fω=fRω+fIωi=0.The real and imaginary terms in expansion equation (25) are arranged into the real part equation fRω and imaginary part equation fIω:(26)fRω=R4ω4+R3ω3+R2ω2+R1ω+R0=0,(27)fIω=I3ω3+I2ω2+I1ω+I0=0.These two equations are quartic equations, and their roots have to beω≠0. Therefore, the imaginary part equation has no constant term and can be written as cubic equation (27).The coefficients of real and imaginary part equations are written as follows:(28a)R4=Imh+ρIB2H4∗+ρmhB4A3∗+ρ2B6A3∗H4∗−ρ2B6A4∗H3∗+ρ2B6A1∗H2∗−ρ2B6A2∗H1∗,(28b)R3=2ρB2IξαωαH1∗+2ρB4mhξhωhA2∗,(28c)R2=−Ikheq−mhkαeq−ρB4kheqA3∗−ρB2kαeqH4∗−4Imhξhξαωhωα,(28d)R1=0,(28e)R0=kheqkαeq,(28f)I3=ρIB2H1∗+ρmhB4A2∗−ρ2B6A1∗H3∗+ρ2B6A3∗H1∗+ρ2B6A2∗H4∗−ρ2B6A4∗H2∗,(28g)I2=−2mhIωhξh−2mhIωαξα−2ρIB4ωαξαH4∗−2ρmhB4ωhξhA3∗,(28h)I1=ρB4kheqA2∗−ρB2kαeqH1∗,(28i)I0=2kheqIωαξα+2kαeqmhωhξh.In this section, the method to solve equivalent linearized equations (17) and (18) is discussed in which the equivalent stiffness coefficients kheq and kαeq are assumed as constant temporarily.In coefficient matrix (24), it contains the unknown ω and eight flutter derivatives Hi∗ and Ai∗i=1,2,3,and4. Since flutter derivatives Hi∗ and Ai∗ are functions of reduced wind speed Vr, expanding equation (25) or real and imaginary part equations (26) and (27) can be regarded as functions of the independent variable ω which includes parameter Vr.In order to find the solution of equations (26) and (27), the following calculation process is required.To do a cycle solving progress, first, a reduced wind speed value is set to obtain eight flutter derivatives and make the real and imaginary part equations onlyω unknown. Then, the positive roots ωR and ωI can be obtained by solving the two equations, respectively. By increasing the reduced wind speed value from zero gradually, several sets of solutions ωR and ωI related to the reduced wind speed can be obtained. Next, the solutions are drawn on the grid in Figure 4; the value of reduced wind speed as abscissa and the value of ωR and ωI as ordinate. Then, all the dots of ωR and ωI are connected to be two lines, and they intersected at a point, which means that the real and imaginary part equations have a same solution; the frequency is ωC=ωR=ωI, and the reduced wind speed is VCr=Vr. So, when Vr=VCr, ω=ωC is the solution of expanding equation (25). Finally, Vcr and ωc are substituted into the coefficient matrix of homogeneous equations (21) and (22); the ratio of vertical to torsional amplitude, H/A, is obtained. In this case, the coefficient determinant equals to zero, and the solutions H and A of homogeneous equations are not equal to zero. It was the solution of the equivalent linearized flutter equations.Figure 4
The intersection of the real part root and imaginary part root.About the intersection point, the value of abscissa is the critical flutter reduced wind speedVCr=10.088, and the value of ordinate is the critical flutter vibration frequency ωC=0.178. Then, the critical flutter wind speed has been calculated that U=75.114 according to Vr=U/Bω, and the vertical-torsional ratio is H/A=9.443, that is, the critical state of the flutter system in this paper. In other words, when the wind speed is 75.114m/s, the suspension bridge will flutter and be instable to bifurcate, and its frequency of vibration is ω=0.178, and its torsional-vertical ratio is H/A=9.443.


## 3.3. Solution of Equivalent Linearization Equations

In this section, the equivalent stiffness coefficientskheq and kαeq in equations (15) and (16) will be taken into account in the solution process because the coefficient kαeq is the function of torsional amplitude A.First, a torsional amplitudeA is the limit cycle amplitude of the system when the flutter occurs is assumed, so the equivalent stiffness coefficient kαeq is determined; the following calculation steps are the same as in Section3.2.The critical flutter solutions including reduced wind speedVCr, frequency ωC, wind speed UC, vertical-torsional ration H/A, and vertical amplitude H can be all obtained. So, these solutions can be interpreted as that the nonlinear flutter system will cause limit cycle oscillation (LCO) when the wind speed is UC. The vertical and torsional amplitude of limit cycle oscillation are H and A, respectively, and the frequency of vibration is ωC.Based on the above analysis and calculation, the characteristic curves of the nonlinear flutter system are drawn in Figures5–8. In order to validate the equivalent linearization method (ELM), numerical solutions by using the Runge–Kutta method [46, 47] are obtained, whose results are shown in Figures 5–8 too.Figure 5
Curve of wind speed and amplitudes of LCOs.Figure 6
Curve of wind speed and frequency of LCOs.Figure 7
Curve of frequency and amplitudes of LCOs.Figure 8
Curve of reduced wind speed and vertical-torsional ratio.As shown in Figure5, the values of wind speed U are taken as abscissa, and the values of vertical and torsional amplitudes are taken as ordinate; then, the same, as in the abovementioned paragraph, two curves in the grid can be obtained. The nonlinear system has flutter instability when critical flutter wind speed is UC=75.114m/s and system bifurcates to LCOs at the same time. Beyond the critical flutter wind speed, the amplitude of LCOs in vertical and torsional degrees increases monotonously when the wind speed is increasing. Also, the solution by ELM agrees well with the numerical results by the Runge–Kutta method.As shown in Figure6, the values of wind speed U are taken as abscissa, and the values of flutter LCOs’ frequency ω are taken as ordinate, and the nonlinear system reaches the critical state of flutter at wind speed UC=75.114m/s, and the critical frequency of flutter is ωC=0.1786. Then, the limit cycle frequency increases monotonically with the wind speed after the bifurcation. Compared with the numerical solution of the Runge–Kutta method at the “×” point in the figure, the frequency of the limit cycle oscillation obtained by the ELM is the larger one.As shown in Figure7, the values of flutter LCOs’ frequency ω are taken as abscissa, and the values of vertical and torsional amplitudes are taken as ordinate; then, two curves are lined through dots of the vertical and torsional degree, respectively, in the grid. The frequency of the LCO of the nonlinear system at the bifurcation point is the critical frequency ωC=0.178. Beyond the critical flutter value, the amplitude of LCOs in vertical and torsional degrees increases monotonously with increasing frequency of LCOs. And, it is that the numerical results are slightly larger than the solutions by ELM.As shown in Figure8, the values of reduced wind speed are taken as abscissa, and the values of the vertical-torsional ratio are taken as ordinate; one curve is lined through the solid square dots. The nonlinear system has flutter instability when critical flutter reduced wind speed VCr=10.088. The vertical-torsional ratio of LCOs decreases monotonously with reduced wind speed.Finally, to validate the feasibility of ELM, the vertical amplitude, torsional amplitudes, and frequency, which at the reduced wind speedVr=10.248, are substituted into the assumed solution h=Heiωt and α=Aeiωt in Section 3.3. Now, the LCO can be drawn into the phaser diagrams shown in Figure 9. The abscissa is the displacement of vertical and torsional degrees, and the ordinate is velocity of vertical and torsional degrees. The solid lines are solutions of ELM, and dots are numerical results. The comparisons show that the results of ELM and numerical method are very closed in the displacement phase, which indicated that the first-order harmonic is dominant in LCO and that the ELM considering the first-order harmonic is practicable. The LCO, which is assumed to be harmonic vibration, is very close to the numerical results. Compared with the LCO of the vertical and torsional degree, the LCO of the vertical degree is closer to the numerical results in displacement and velocity, and the velocity in LCO of the torsional degree is larger than the numerical results with an error of 16% which may be related to the removal of higher-order harmonic vibration in ELM.Figure 9
Limit cycle oscillation (LCO) vs. numerical solutions. (a) LCO in the vertical degree. (b) LCO in the torsional degree.
(a)(b)Therefore, the flutter bifurcation problem of the nonlinear system can be analyzed by the calculation result of the ELM.


## 4. Analysis of the Hopf Bifurcation

Hopf bifurcation refers to the phenomenon that, as the system parameters change to exceed the critical value, the system suddenly moves flutter and periodic oscillation from a stable equilibrium state before. From the phase diagram, when beyond the critical flutter value, the stable focus of the system will suddenly appear as an isolated closed trajectory called limit cycle oscillation (LCO).In the solution of the previous section, when the determinant of the flutter system is zero, the homogeneous equation has nonzero solutions. And, the flutter system is in a critical state, that is, when the system is slightly disturbed or exceeds the critical state, there will be a bifurcation and a nonzero solution will appear at this time.In the vicinity of the flutter critical value, two reduced wind speed values areVr=10.0866 and Vr=10.0966, and the flutter critical frequency ωC=0.1786 is selected; then, the LCOs are solved by using the Runge–Kutta method and are drawn on the phase diagram.As shown in Figures10(a) and 10(b), when the reduced wind speed is less than the critical value that Vr=10.0866, the trajectory of the phase diagram has an amplitude of 5×10−4, which can be regarded as a stable focus, and the system does not vibrate. When the reduced wind speed Vr=10.0906 that is higher than the critical value, results in Figures 10(c) and 10(d) show that vertical amplitude of LCO is 0.035 and that of the torsional amplitude is 0.0037, which is a stable limited amplitude vibration, forming an isolated closed trajectory. When the parameter Vr gradually increases, the nonlinear system bifurcates from the stable state to the LCO state. So, the bifurcation of the system in the critical state belongs to the supercritical Hopf bifurcation.Figure 10
LCOs of two different reduced wind speed values. (a) LCO whenVr=10.0866. (b) LCO when Vr=10.0866. (c) LCO when Vr=10.0906. (d) LCO when Vr=10.0906.
(a)(b)(c)(d)


## 5. Effect of Damping Ratio and Cubic Stiffness on Bifurcation and LCOs

It is known from the previous section that a suspension bridge system with strong nonlinearity undergoes a supercritical Hopf bifurcation in the flutter critical state and has LCOs appeared on the right side, while the damping and nonlinear stiffness coefficients of the system may affect the bifurcation and LCOs. In this section, the effect of damping ratio and cubic stiffness on bifurcation and LCOs of the nonlinear system is taken into discussion. Based on the original bridge nonlinear flutter system, sets of damping ratios for vertical and torsional degrees and different cubic stiffness coefficients were set, respectively, to discuss the effects on bifurcation state and LCOs.The damping ratio of a steel integral structure generally ranges from 0.003 to 0.007, and it will change with service life and different construction conditions. Also, tuned mass dampers are usually installed to adjust the bridge damping ratio [48, 49]. Therefore, by analyzing the effects on critical flutter state and LCOs when the bridge system damping ratio is changed, we can get the degree of freedom damping ratio change that has the greatest influence on bridge flutter.Three working conditions are discussed, as shown in Table3. Condition 1 is the damping ratio of the original bridge system. In condition 2, the torsional damping ratio is increased by 40% to 0.007 on the basis of condition 1, and the vertical damping ratio remains unchanged, and the effect of changing the torsional damping ratio is studied. In condition 3, on the basis of case 2, the vertical damping ratio is reduced by 40% to 0.003, and the torsional damping ratio remains unchanged. The effect of the vertical damping ratio is studied.Table 3
Values of the damping ratio in three conditions.
Condition 1Condition 2Condition 3Vertical damping ratioξh0.0050.0050.003Torsional damping ratioξα0.0050.0070.007Using the ELM and solving the equivalent linearized equations in this paper, the bifurcation and LCOs can be calculated, and the results are shown as follows.As shown in Figure11, the nonlinear systems of the three conditions reach the critical state of flutter at wind speed U=75.1144, U=76.5683, and U=76.5602, and bifurcation occurs. LCOs appear on vertical and torsion degrees simultaneously. Conditions 2 and 3 bifurcate later than condition 1, and their bifurcation curves almost coincide. The amplitude of the limit cycle of the three conditions changes with the increase of wind speed similarly. For the vertical degree, the amplitude of the LCO in the higher reduced wind speed of the bifurcation of condition 2 and condition 3 is slightly larger than that of condition 1. As for the torsional degree, the amplitude of the limit cycle is always smaller than that of condition 1 when conditions 2 and 3 have the same wind speed.Figure 11
Curves of wind speed and amplitudes of LCOs in three conditions. (a)U vs. H. (b) U vs. A.
(a)(b)In Figure12, all three conditions were flutter bifurcation and their flutter critical frequency were ω=0.1768, ω=0.1773, and ω=0.1773. When beyond the flutter critical state, their amplitudes of LCOs were increasing by the growing of wind speed and frequency. In both vertical and torsional degrees, the amplitudes of LCOs in first condition are larger than the other two.Figure 12
Curves of frequency and amplitudes of LCOs in three conditions. (a)ω vs. H. (b) ω vs. A.
(a)(b)As shown in Figure13, at the moment of the bifurcation of the critical flutter state, the critical wind speed of conditions 2 and 3 is larger than condition 1. In addition, the critical flutter wind speed of condition 2 is slightly greater than that of condition 3, and the critical flutter frequency is smaller than condition 3. The wind speed and frequency relationship lines of conditions 2 and 3 almost coincide with each other. The frequency of LCO in all three conditions shows a linear increasing trend with the increases of wind speed.Figure 13
Curves of wind speed and frequency of LCOs in three conditions.As shown in Figure14, among the three working conditions, after the bifurcation of the nonlinear system, the ratio of the amplitude of the LCO of the vertical to the torsional decreases as the reduced wind speed increases. The vertical-torsional ratio of condition 2 and working condition 3 is greater than that of working condition 1 and shows a trend of rapid decrease.Figure 14
Curves of reduced wind speed and vertical-torsional ratio in three conditions.According to the above analysis, comparing conditions 1 and 2, which only changing the torsional damping ratio, have an obvious influence on the Hopf bifurcation point and the LCOs after bifurcation of the nonlinear system. Increasing the torsional damping ratio, the Hopf bifurcation point of the nonlinear system, that is, the delay of the critical state, bifurcation occurs at a larger reduced wind speed and real wind speed. In the bifurcation critical state, the critical flutter frequency is smaller, and the amplitude and frequency of the LCOs are smaller after the bifurcation, and the vertical-torsional ratio of the LCOs amplitude is larger. Comparing conditions 2 and 3, which only changing the vertical damping ratio, have little effect on the Hopf bifurcation point and various parameters after bifurcation of the nonlinear system. If the torsional damping ratio of the suspension bridge decreases, the flutter critical wind speed will seriously reduce, but also the LCO amplitude of flutter will grow up.Then, we discuss another set of comparisons, the influence of different cubic stiffness coefficients on the flutter bifurcation and limit cycles of the system. Assuming that the nonlinear stiffness coefficient of the bridge vibration system is set to 1.5 times the initial stiffness coefficient, that is,eα=1500mhωα2, and called the “strong stiffness” system. The damping ratio of the vertical and torsional degree remains unchanged, that is, ξh=0.005 and ξα=0.005; the calculation result is shown as the following figures. The dotted line is the initial system of the original system, and the solid line represents the “strong stiffness” system with 1.5 times the cubic stiffness coefficient.As shown in Figures15 and 16, both systems reach a critical state at the same time and flutter occurs. When the critical reduced wind speed VCr=10.088 and the critical wind speed UCr=75.114m/s, the system is in a critical state, the critical vibration frequency is ωC=0.178Hz, and the limit cycle vibration appears. The amplitude of the limit cycle oscillation increases with the increase of reduced wind speed and real wind speed. As the amplitude increases, the vibration frequency also increases. The “strong stiffness” system is 1.5 times the cubic stiffness coefficient at the same reduced wind speed, wind speed, and frequency. The limit cycle oscillation amplitude is smaller than the initial basic system. It can be seen from Figures 15 and 16 that the “strong stiffness” system is at a higher wind speed and the vibration frequency is higher when the amplitude of the limit cycle of the two systems is same.Figure 15
Curve of wind speed and amplitudes of LCOs.Figure 16
Curve of frequency and amplitudes of LCOs.As shown in Figures17 and 18, in the relationship between the wind speed vs. frequency of the limit cycle and the reduced wind speed vs. vertical-torsion ratio, the curves of the two systems coincide on the same line.Figure 17
Curve of wind speed and frequency of LCOs.Figure 18
Curve of reduced wind speed and vertical-torsional ratio.According to the above discussion of the two systems with different cubic stiffness coefficients, they both reach the critical state of flutter at the same wind speed value. The larger cubic stiffness coefficient has a positive aspect to the limit cycle vibration of the flutter system. In the same wind speed, its limit cycle oscillation amplitude is smaller. When the limit cycle vibration of the same amplitude occurs, the vibration frequency is higher. The relationship between limit cycle oscillation frequency and wind speed and vertical-torsion ratio and reduced wind speed has not changed.


## 6. Results and Discussion

In this paper, a bridge flutter model with cubic nonlinear stiffness was taken to study, which reflects the subject that the large-span characteristic of the Nizhou Waterway Bridge. The qualitative and quantitative flutter analysis for a long-span suspension bridge with strong nonlinearity in structural stiffness is carried out by an approximate analytical equivalent linearization method, revealing the limit cycle oscillation phenomena of bifurcation beyond the flutter critical state.When the nonlinear flutter system reaches a critical state, a supercritical Hopf bifurcation occurs and LCOs are born to the right of the bifurcation point. This nonlinear flutter system will vibrate periodically with a certain amplitude and frequency at a certain wind speed beyond the critical speed. The amplitude and frequency of the vibration are increasing as wind speed grows up.Also, analyzing the influence of structural factors on them, this provided a reference basis for wind resistance design of long-span suspension bridges. In particular, increasing the torsional damping ratioξh had effectively improven the critical flutter wind speed of the bridges. On the contrary, the decrease of the torsional damping ratio in a service bridge caused by structural factors change will seriously reduce the critical flutter wind speed. However, increasing the vertical damping ratio ξh had little effect on the delay of the critical flutter state, the amplitude and frequency of LCOs after bifurcation, and the vertical-to-torsion ratio, which is not so apparent as the change is brought about by the increase of the torsional damping ratio. At the same time, the use of structural measures to improve the torsional stiffness of the structure also had a positive effect. A larger cubic stiffness coefficient took a positive effect on the bridge flutter system although the critical state had no difference, but the vibration amplitude of the limit cycle oscillations at one wind speed was reduced. When the limit cycle oscillation occurs with the same amplitude, a system with a higher cubic stiffness coefficient had a higher vibration frequency.



---

*Source: 5572429-2021-07-20.xml*