On boundary damping to reduce the rain–wind oscillations of an inclined cable with small bending stiffness

In this paper, a model will be derived to describe the rain–wind-induced oscillations of an inclined cable. Water rivulets running along the cable and aerodynamics forces acting on the cable are taken into account to describe these oscillations. A boundary damper is assumed to be present near the lower endpoint of the cable. For a linearly formulated initial-boundary value problem for a tensioned beam equation describing the in-plane transversal oscillations of the cable, the effectiveness of this damper is determined by using a two-timescales perturbation method. It is shown how mode interactions play an important role in the dynamic behaviour of the cable system. Some resonant and non-resonant cases have been studied in detail.


Introduction
The study on how to damp vibrations in stay cables, which are attached to a pylon at one end and to a bridge deck at the other end, is of great importance not only in structural engineering but also in applied mathematics. The combined effect of rain and of wind can change the aerodynamic properties of the cablestayed bridge and can lead to relatively large amplitude vibrations of the cables. For example, one can refer to the Erasmus bridge in Rotterdam, which started to vibrate heavily under mild wind-rain conditions shortly after its opening in 1996. This bridge was temporarily closed to the traffic as a safety precaution. As a temporarily measure, polypropylene ropes were installed between the cables and the bridge desk. Later, these ropes were replaced by hydraulic dampers as a permanent measure and by conforming to the aesthetic of the bridge designed by the architects. Much research on this problem has been done both numerically and experimentally to understand the mechanisms of rainwind-induced vibration of inclined stay cables [10,11]. For more recent experiments and numerical computations, the reader is referred to the following articles [5,8,13,18].
As has been observed from wind-tunnel experiments, raindrops hitting the inclined stay cable cause the generation of one or more rivulets on the surface of the cable. The presence of flowing water on the cable changes the mass and the aerodynamic properties of the bridge system that can lead to instabilities. These water rivulets on the cable surface can be considered as a time-varying mass in the system [1,3]. In order to suppress undesirable vibrations in bridge structures, different kinds of dampers such as tuned mass dampers and oil dampers can be installed between the cables and the bridge deck. Numerous work has been done from the theoretical and the experimental point of view to predict the optimal damper location and type. Jacquot [15] calculated the optimal value and location of the viscous damper which is located in a randomly forced horizontal cantilever beam.
These physical problems of rain-wind-induced vibration of inclined cables can be modelled mathematically by initial-boundary value problems for string-like or beam-like problems. For string-like problems, the static state due to gravity and the dynamic state due to a parametrical and a transversal excitation at one of the ends of the inclined string were studied in [7]. For the inclined cable subjected to wind with a moving rivulet on its surface, the nonlinear dynamic model is investigated by considering the equilibrium position of the rivulet [19]. The effect of the static condensation of the longitudinal displacement due to the cable inclina-tion and the cable total tension is investigated by using numerical and analytical techniques in [25]. In addition, the interaction among the three excitation sources: self excitation, which is caused by a mean wind flow, and external and parametric excitations due to vertical motion of the ground support, on the nonlinear dynamics of the inclined cable related to a cable-stayed bridge were studied in [20]. In many papers, the rain-windinduced vibrations of inclined taut cables have been studied, but much less work has been done on the rainwind-induced vibrations of an inclined cable subjected to wind with time-varying rivulets on its surface.
From the practical point of view, the viscoelastic materials are modelled with a combination of springlike and dashpot-like elements at the boundary in order to represent restorative force component and damping component, respectively. In this paper, we consider the linear (and nonlinear) dynamic response due to only a damper near one of the support ends of the inclined cable. On the cable, a time-varying (rain) mass rivulet is assumed to be present (see also Fig. 1). The stationary wind flow will lead to the phenomenon of self excitation of the cable. The damper location in the cable-stayed bridge is quite close to the anchorage of the cable. For more information on the effect of the bending stiffness of a tensioned cable for varying damper locations, the reader is referred to [14,21,22].
The aim of this paper is to provide an understanding of how effective boundary damping is for inclined stretched beams with a small bending rigidity. These problems for strings or beams are considered to be basic models for oscillations of cables from a practical viewpoint. The outline of this paper is as follows. In Sect. 2, we apply a variational method in order to derive the governing equations of motion of the tensioned Euler-Bernoulli beam, and obtain a system of three coupled PDEs. By using Kirchhoff's approach, the number of PDEs in the system is reduced to two PDEs (one for the in-plane motion, and one for out-ofplane motion) or to a single PDE (only for the in-plane motion). The main aim of this section is to give a model to describe the dynamics of rain-wind-induced vibrations of an inclined beam, where the gravity effect is schematically shown in Fig. 1. In Sect. 3, we only consider the in-plane motion of the inclined beam with only one rivulet on the cable. In Sect. 4, the two-timescales perturbation method is used to solve the problem and some (non)resonance frequencies are determined. In this paper, we will consider the pure resonance case and the non-resonance case. All other resonance cases can be investigated similarly. Finally, the conclusions are presented in Sect. 5.

Equations of motion
We consider an inclined, perfectly flexible, elastic beam with a small bending rigidity on a finite interval x ∈ [0, L], which is attached to a dashpot λ at x = L, and assumed to be simply supported at x = 0 (see Fig. 1). u, w and v are the displacements in x-direction, ydirection and z-direction, respectively. From Newton's second law, the equations of motion can be stated as follows: the time derivative of the linear momentum of the system is equal to the sum of external forces which are elastic forces, gravity, drag and lift forces due to the uniform wind flow v 0 , blowing under a yaw angle β as can be seen in Fig. 1. It is assumed that the water rivulet is not blown off the cable. Furthermore, we assume that the tension due to stretching in the beam is large enough, such that the small sag of the beam, that is, the displacements in x and z-directions due to gravity, can be neglected. Similarly, due to this tension, we neglect the wind force along the cable in x-direction. Hence, the total tension in the beam at x = x 0 can be written by: where T 0 is the pretension in the beam, g is the acceleration due to gravity, A is the cross-sectional area of the beam, α is the angle between the beam and the horizontal plane, and M(x, t)) is the mass density of the beam including the time-varying water rivulet. Let the coordinates (x, y, z) of a material point of the unstretched beam be (x, 0, 0) with x ∈ [0, L], where the x-axis, y-axis and z-axis are defined in Fig. 1. The dynamic displacement of this material point is denoted by u(x, t)i, v(x, t)k and w(x, t)j, where i, k and j are the unit vectors along the x-axis, y-axis and z-axis [7]. The vector position R(x, t) of this material point in the dynamic state is obtained by: where E is Young's modulus. The relative strain per unit length of the stretched beam is written by: where u x , v x ,and w x represent the derivative with respect to x of u(x, t), v(x, t) ,and w(x, t), respectively. By assuming that u 2 x is small with respect to u x and v 2 x + w 2 x , and by expanding the square-root in a Taylor series, we have approximately For the curvature of the beam in the (x, y, z)-space (out of plane) or for the curvature of the beam in the (x, y)-plane (in plane), the reader is referred to [6,23]. The axial strain of a generic point of the beam located at distances y and z (in the y and z directions) from the centre line of the beam is defined by: Hence, the total strain of a line-element of the beam is approximately given by: The equations of motion describing the vibrations of the beam will be obtained by using Hamilton's principle [12]. In order to apply this principle, the kinetic energy and potential energy for the beam should be defined and calculated. The potential energy of the beam is given approximately by (using Hooke's Law [16]) and defined as: (7) where I y and I z represent the axial moments of area about the y-and z-axes, respectively. In addition, the kinetic energy of the beam is given by: The , then the three Euler-Lagrange equations of the variational problem δF = 0 , which we have to solve according to the Hamiltonian principle, are as follows: or equivalently, the equations of motion are given by where F y and F z are the aerodynamic forces in the in-plane and in the out-of-plane, respectively, and are defined by Here D and L are the magnitudes of the drag and lift forces, respectively, which may be given by: where ρ a is the air density, d is the diameter of the cross section of the circular part of the beam cable, L is the length of beam, v s is the virtual wind velocity which is given by  Figure 2 shows the centre of a cross section of the cable with rivulets. ψ 1 and ψ 2 are the displacements of the upper and the lower rivulet, respectively, and φ * i is the wind attack angle, which is defined by φ * i = φ −ψ i , i = 1, 2. The angle between the virtual wind velocity v s and the horizontal axis is defined by φ, which can be expressed as: where v ∞ = v 0 cos 2 (β) + sin 2 (α)sin 2 (β) is the effective wind flow, and γ is the angle of attack, which depends on the inclination angle α and the yaw angle β defined by Geurts and van Staalduinen [11] γ (α, β) = arcsin sin(α)sin(β) As can be seen in Fig. 2, the position of the upper and lower rivulet is determined by ψ 1 and ψ 2 , respectively. From experimental data [4], we may assume that the total time-varying mass due to these rivulets on the cable changes periodically and may be defined to have the form [1,2] where M 1 > 0 and M 2 > 0 are the constant mass of the upper and lower rivulets, respectively. A 1 > 0 and A 2 > 0 are small parameters, and γ 1 > 0, γ 2 > 0, If M 0 > 0 is the constant mass of the beam cable, the total mass can be given by: where i ) coefficients may be obtained from wind-tunnel measurements. A quasi-steady assumption can be used due to the low-frequency oscillations of the cable. Yamaguchi shows the experimental results for the drag, lift and moment coefficients at various angles of wind attack with the ratio of diameters of rivulet and cable in [27]. These results show that the drag, lift and moment coefficients for different ratios of diameters of rivulet and cable are similar to each other. Due to these results, we define that the drag and lift coefficients may be written in periodic form in a small α 0 , β 0 and γ 0 neighbourhood of fixed "rivulet profiles as follows: Here κ i ,κ i andκ i are constants for i = 1, 2 indicating the effect of the mass change of the rivulets on the drag and lift coefficients [17], and where α 0 , β 0 and σ 0 are small parameters, and γ 0 > 0, Ω 0 > 0. By neglecting terms of degree four and higher terms, we obtain the in-plane wind force as follows: and similarly, we compute the out-of-plane wind force as follows: where A i j and B i j for i, j = 0, 1, 2, 3 as follows: and The detailed expressions of the a i jk -and b i jkcoefficients for k = 1, 2, 3, 4 can be found in "Appendix A". When we assume only in-plane horizontal cable motion with only the upper rivulet present, that is, γ (α, β) = 0, μ 2 (x, t) = 0, and κ 1 =κ 1 =κ 1 = μ 1 (x, t) = 1, the same coefficients of c D , c L 1 and c L 3 are obtained as in [26].

Further simplifications
For further calculations, we only consider in-plane motion (w = 0) of the inclined beam, and we assume that there is only one rivulet on the cable, that is, whereÃ is a small parameter. We substitute Eqs. (27) and (28), into Eqs. (10)- (11) and rewrite in order to obtain the equations of motion for the beam with a time-varying mass, yielding where A i j for i, j = 0, 1, 2, 3 are defined in Eq. (29). As can be seen in Fig. 1, the inclined cable is attached to a sliding damper at x = L and to a pylon at x = 0. The boundary conditions are given by: and and the initial conditions are Equations (31)-(32) represent the in-plane motion of the inclined beam cable system in the longitudinal and transversal direction, that is, in x-, and y-direction. We introduce the new variablesū(x, t) andv(x, t) as defined by: whereû(x) andv(x) are the "stationary solutions by neglecting small-time perturbation (see "Appendix B"). When we substitute the new variables Eqs. (37)-(38) into the Eqs. (31)-(32), we obtain These coupled partial differential equations can be reduced to a single partial differential equation by applying Kirchhoff's approximation. It will be assumed where is a small parameter with 0 < 1. Then, by using these assumptions, Eq. (39) up to order becomes and then from 0 to L, obtaining Hence, When we substitute Eq. (44) into Eq. (42), we obtain Similarly, the equation for v in Eq. (40) can be rewritten in where from Eq. (45) and û x +v 2 x 2 from Eq. (241) as given in "Appendix B" into Eq. (46) we obtain In order to put the equations in a non-dimensional form, the following dimensionless quantities are used: Then, Eq. (47) in a non-dimensional form becomes with the boundary conditions and the initial conditions We also assume thatv(x, t) = v(x, t),η 10 = η 10 ,η 2 = η 2 , η 3 = η 3 , andÃ = σ . The asterisks indicating the dimensional quantities are omitted in Eqs. (48) through (51), and henceforth for convenience. Then, by using these assumptions, Eq. (48) up to order 2 becomes with the boundary conditions and the initial conditions In the following section, the initial-boundary value problem Eqs. (52)-(55) will be studied further.

Application of the two-timescales perturbation method
In this section, the initial-boundary value problem Eqs. (52)-(55) will be studied and an approximation of the solution of the initial-boundary value problem up to order will be constructed by using a two-timescales perturbation method. We assume that v(x, t) can be expanded in a formal power series in , that is, where all v i (x, t) and their derivatives for i = 0, 1, 2, · · · are O(1) on a time scale of order −1 . The approximation of the solution may have secular terms which are unbounded terms in time. In order to avoid these secular terms, we will apply the two-timescales perturbation method by introducing a slow time scale The following transformations are needed for the time derivatives v t = y t + y τ , (58) v tt = y tt + 2 y tτ + 2 y τ τ . (59) Substitution of Eqs. (57)-(59) into Eq. (52) yields: with the boundary conditions and the initial conditions Assuming that then by collecting terms of equal powers in , it follows from the problem for y(x, t, τ ; ) that the O(1)problem is and that the O( )-problem is The method of separation of variables will be applied to solve the problem Eqs. (65)-(69). The solution of the O(1)-problem may be given in a special form By substitution of Eq. (75) into Eq. (65) and by dividing the so-obtained equation by T (t, τ )φ(x) yields: A separation constant is defined −ω so that the timedependent part of the product solution oscillates for real and positive eigenvalues (for the proof, we refer the reader to [24]). We obtain a time-dependent part and the general solution of the time-dependent part is a linear combination of sines and cosines in t, where σ 1 and σ 2 are arbitrary function in τ . In addition, we obtain a space-dependent part and the boundary conditions Eqs. (66)-(68) yield The characteristic equation for Eq. (79) is given by . The non-trivial solutions are found by using the boundary conditions Eq. (80), leading to the characteristic equation It follows from Eq. (83) that the eigenvalues ω n = μb 4 n +b 2 n can be numerically computed for given values of μ. The first ten eigenvalues ω n are listed in Table 1. The eigenfunctions belonging to different eigenvalues are orthogonal with respect to the inner product; for details, the reader is referred to [24], where θ n = − b n cos(b n ) a n cosh(a n ) , a n = 1+ √ 1+4μω n 2μ , and b n = −1+ √ 1+4μω n 2μ . Hence, infinitely many non-trivial solutions of the initial-boundary problem Eqs. (65)-(69) have been determined. By using the superposition principle, the solution is obtained where φ n (x) is given by Eq. (85), and where A n and B n are arbitrary functions in τ which can be used to avoid secular terms in y 1 (x, t, τ ). By using the initial conditions Eq. (69), we can determine A n (0) and B n (0), which are given by: and where Now, we solve the O( )-problem Eqs. (70)-(74). Due to having an inhomogeneous boundary condition Eq. (73), we use the following transformation to convert the problem into a problem with homogeneous boundary conditions where Substituting Eq. (90) into Eqs. (70)-(74), we obtain where h(t, τ ) is given by Eq. (91), and where h(0, 0), and h t (0, 0) are given by In order to solve Eqs. (92)-(96), V (x, t, τ ) is written in the following eigenfunction expansion and by substituting Eq. (99) into the partial differential equation Eq. (92), we obtain We expand x 4 −2x 3 12μ+2 and 12x 2 −12x−24μ 12μ+2 into a series of eigenfunctions φ m (x), and we obtain where and where ζ m is given by Eq. (89). By multiplying both sides of Eq. (100) by φ n (x), and then by integrating from x = 0 to x = 1, and by using the orthogonality properties of the eigenfunctions, we obtain V n tt (t, τ ) + ω n V n (t, τ ) = −c n h tt (t, τ ) + d n h(t, τ ) − 2T n tτ (t, τ ) where and

It follows from Eqs. (91) and (86) that h(t, τ ) and h tt (t, τ ) can be written as
Hence, by using Eqs. (113) and (114) it follows that Eq. (105) can be rewritten as 10 √ ω n +λφ n (1) √ ω n (c n ω n + d n ) It can be easily seen from Eq. (115) that there are infinitely many values of Ω 1 , which can cause internal resonances. The possible resonance cases are as follows: (i) The non-resonant case: if Ω 1 is not within an order -neighbourhood of the frequencies √ ω n ± √ ω m for all m and n.
(ii) The near resonance case (Resonance detuning): if Ω 1 is within an order -neighbourhood of √ ω n ± √ ω m for certain fixed m and n.
(iii) The pure resonance case: if Ω 1 = √ ω n + √ ω m for certain fixed m and n. This is the sum-type resonance case. If Ω 1 = √ ω n − √ ω m for certain fixed m and n, it is the difference-type resonance case. In the following subsections, we will study the non-resonant case, the sum-type resonance case, the difference-type resonance case and by considering only some of the first few modes and omitting the higherorder modes.
For simplicity, we will now assume that . Resonance can occur when Ω 1 = ± √ ω n ± √ ω m for some n and m. We have to consider the following three cases: which has no solution since the right-hand side is negative, while the left-hand side is positive.
, which has only solutions by looking at Table  4.1 for certain values of μ. There are possibilities that three modes are interacting. For example,

The non-resonant case
When Ω 1 = √ ω n ± √ ω m + O( ) for all m and n only resonances occur due to the fourth and fifth term in the right-hand side of Eq. (115), that is, Eq. (115) can be rewritten as: 10 √ ω n +λφ n (1) √ ω n (c n ω n + d n ) where "NST" stands for terms that lead to nonsecular terms in V n . In order to remove secular terms, it follows from Eq. (116) that A n (τ ) and B n (τ ) have to satisfy where X n and Y n are defined by The solution of Eqs. (117) and (118) is given by where A n (0) and B n (0) are given by Eqs. (87) and (88), respectively. By using Eqs. (103) and (104) with Eqs. (79) and (80) it can easily be shown that 1 2 (c n ω n + d n ) = μ ζ n (12μ+2) (1). Hence, Eq. (106) with Eq. (119) and Eq. (122) can be rewritten as: Now, from Eq. (116) with Eqs. (117)-(123), we can obtain V n (t, τ ) straightforwardly. Obviously, V n (t, τ ) will be bounded on a time scale of order 1/ , and so will be V n (x, t, τ ).
In Table 2, numerical approximations of ω n and the damping parameterp n are given for different values of μ. For different bending stiffness, we can choose the damping parameter in such a way that all modes are damped by takingλ appropriately compared to the η 10 coefficient. It is also clear to see from Table 2 that for smaller values of μ, we should take largerλ to have damping for all modes. In Table 3, it can be seen as expected that the bending stiffness μ and the damping parameterλ influence the stability.
For μ = 0.001, Figs. 3 and 4 demonstrate the vibration response at the middle x = 0.5 and at the end x = 1    Fig. 3 is plotted forλ = 0.5 and Fig. 4 is plotted forλ = 2; the initial conditions are specified as v 0 (x) = 0.1sin(π x) and v 1 (x) = 0.05sin(π x). For μ = 1, Figs. 5 and 6 illustrate similar behaviour. These figures show that the amplitudes of the transverse vibrations decrease faster for increasing of μ andλ. Similar results for the string-like problem have been observed in [9].

4.2
The sum-type resonance case: We will consider now only the first two modes and omit the higher-order modes. That is, we will assume that for the given μ values only these two modes might occur in a resonance interaction. Equation (115) can be rewritten as and a similar equation for V 2 (t, τ ) can be obtained from Eq. (124). In order to avoid secular terms, it follows from the equations for V 1 (t, τ ) and V 2 (t, τ ) that A 1 (τ ), B 1 (τ ), A 2 (τ ) and B 2 (τ ) have to satisfy Table 3 Numerical approximations ofω n in the non-resonance case for η 10 = η 2 = γ 1 = 0.1, σ = 2.8, and different values of μ ω nω jλ = 0λ = 0.05λ = 0.1λ = 0.5 12096ĩ 12096i 19285ĩ 26094ĩ 04832ĩ 10884ĩ 15541ĩ where X 1 , X 2 , Y 1 , Y 2 , Z 1 , Z 2 , C 1 and C 2 are defined by 04557ĩ 10610ĩ We obtain from the system Eqs. (125)-(128), and whereẊ represents the derivative of X with respect to τ . The matrix A for a given configuration is only depending on the damping parameterλ. The other parameters are determined by the physics, which are taken from [1,3] in Table 4. In order to stabilise the system, we should determine the damping parameterλ in such a way that all the real parts of the eigenvalues of the matrix A are negative. As can be seen in Table 4, if we assumeλ = 0, there is an instability in the sys-tem due to the wind force η 10 . While for increasing values ofλ, we observe that the cable system becomes stable, which also depends on the value of the bending stiffness μ.

4.3
The difference-type resonance case: A similar analysis, as given for the sum-type resonance case, can also be applied in the difference-type resonance case. We will consider Ω 1 = √ ω 2 − √ ω 1 , which is assumed to have only solution m = 2, n = 1 (or m = 1, n = 2). Then, Eq. (115) can be rewritten as 11951ĩ 04791ĩ 10926ĩ and a similar equation for V 2 (t, τ ) can be obtained from Eq. (129). In order to avoid secular terms, it follows from the equations for V 1 (t, τ ) and V 2 (t, τ ) that A 1 (τ ), where X 1 , X 2 , Y 1 , Y 2 ,Z 1 ,Z 2 ,C 1 andC 2 are defined by: We obtain from the system Eqs. (130)-(133), and whereẊ represents the derivative of X with respect to τ . As in the sum-type resonance case, this matrix A for a given configuration is only depending on the damping parameterλ. In Table 5, it can easily be seen that there is a change from instability to stability, which is aroundλ = 0.05.

Conclusions
In this paper, initial-boundary value problems for a tensioned beam equation are studied. The model is derived to describe the rain-wind-induced oscillations of an inclined cable. We applied a multiple-timescales perturbation method in order to observe whether or not mode interactions between vibration modes occur for certain values of the bending stiffness and the damping parameter. The results show that the system in both the pure resonance case and the non-resonance case can be stabilised by a boundary damper. Some of these cases are studied in Sect. 4. Mode interactions between two and more modes depending on the bending stiffness μ are possible. More complicated resonance cases can be Ω 1 = √ ω n ± √ ω m or Ω 1 = 2 √ ω n , when for some fixed N and M. These cases still have to be studied, and can be studied by using the techniques as shown in Sect. 4 of this paper.
This paper provides an understanding of how effective boundary damping can be for the in-plane transversal oscillations of the cable. The same approach can be used for in-plane and out-of-plane transversal oscillations of elastic structures.