J T Katsikadelis

G C Tsiatas
email: gtsiatas@central.ntua.gr

Nonlinear dynamic stability of damped Beck's column with variable cross-section

Keywords: Beck's column, nonlinear dynamic stability, nonlinear dynamic analysis, analog equation method, variable cross-section

In this paper the nonlinear dynamic stability of Beck's column with variable mass and stiffness properties in the presence of damping (both internal and external) is investigated using a complete nonlinear dynamic analysis. This approach permits the examination of the global stability of the system in contrast to the static nonlinear one, which, though more economical in computational cost, is associated only with the loss of local stability via flutter or divergence. The governing equations describing the dynamic response are derived in terms of the displacements taking also into account the axial deformation, which has a striking influence on the critical load. Since the cross-sectional properties of the beam vary along its axis, the resulting coupled nonlinear differential equations have variable coefficients. Their solution is achieved using the analog equation method (AEM) of Katsikadelis. Besides its accuracy and effectiveness, this method overcomes the shortcoming of a FEM solution, which may experience lack of convergence. Interesting conclusions are drawn. The important, however, finding is that the inclusion of the axial deformation affects highly the critical load of Beck's column with varying cross sectional properties, while it leaves it unaltered for Beck's column with uniform cross section.

Introduction

The influence of damping on the stability of linear elastic nonconservative uniform columns has been extensively investigated. The related main conclusions are the stabilizing effect of the external damping [START_REF] Plaut | The effect of external damping on the stability of Beck's column[END_REF] and to the paradox of the destabilizing effect of the internal damping [START_REF] Bolotin | Nonconservative Problems of the Theory of Elastic Stability[END_REF][START_REF] Herrmann | On the destabilizing effect of damping in nonconservative elastic systems[END_REF]. However, when nonuniform columns are examined [START_REF] Kar | Stability of a nonuniform cantilever subjected to dissipative and nonconservative force[END_REF][START_REF] Rao | Stability of tapered cantilever columns subjected to a tipconcentrated follower force with or without damping[END_REF] the influence of internal damping can have a stabilizing character under certain conditions, while the external damping has always a stabilizing effect [START_REF] Kar | Stability of a nonuniform cantilever subjected to dissipative and nonconservative force[END_REF]. Some recent results for linear viscoelastic models, in which creep, relaxation and hysteresis effects are taken into account, can be found in a survey paper by Gaul [START_REF] Gaul | The influence of damping on waves and vibrations[END_REF]. A special case of the influence of damping on Hamiltonian autonomous systems is very recently reported by Kounadis [START_REF] Kounadis | Hamiltonian weakly damped autonomous systems exhibiting periodic attractors[END_REF].

The stability behavior of nonconservative columns can be examined only via nonlinear analysis, static or dynamic [START_REF] Vitaliani | Finite element solution of the stability problem for nonlinear undamped and damped systems under nonconservative loading[END_REF][START_REF] Kounadis | On the paradox of the destabilizing effect of damping in non-conservative systems[END_REF][START_REF] Kounadis | Non-potential dissipative systems exhibiting periodic attractors in region of divergence[END_REF][START_REF] Bolotin | Effect of damping on the postcritical behaviour of autonomous non-conservative systems[END_REF]. The latter one permits the examination of the global stability of the system instead of the static nonlinear analysis, which is associated with the loss of local stability, via flutter or divergence [START_REF] Kounadis | Some new instability aspects for nonconservative systems under follower loads[END_REF][START_REF] Kounadis | On the failure of static stability analyses of nonconservative systems in regions of divergence instability[END_REF]. The work that has been done on the nonlinear dynamic analysis is limited only to uniform damped Beck's columns. The bifurcation may be subcritical or supercritical. The subcritical can exhibit instability below the critical load (from the linear theory), while the supercritical is the stable solution to which all solutions tend independently of the initial conditions. The initial studies of Kolkka [START_REF] Kolkka | On the non-linear Beck's problem with external damping[END_REF] and Chen [START_REF] Kolkka | On the non-linear Beck's problem with external damping[END_REF] have shown that the bifurcation is supercritical for the uniform Beck's column.

Recently, Andersen and Thomsen [START_REF] Andersen | Post-critical behaviour of Beck's column with a tip mass[END_REF] studying a uniform Beck's column with a tip mass at its free end observed that the rotary inertia of the mass can change the type of bifurcation from supercritical to subcritical. To the authors knowledge publications on the solution of the problem of Beck's column with variable mass and stiffness properties are not available in literature.

A c c e p t e d m a n u s c r i p t 3

In this paper the nonlinear dynamic stability of Beck's column with variable mass and stiffness properties in the presence of damping (both internal and external) using a complete nonlinear dynamic analysis is investigated. For homogeneous material variable mass and stiffness is due to the variation of the beam cross-section. The nonlinearity results from retaining the square of the slope in the strain-displacement relations (intermediate nonlinear theory). In this case the transverse deflection affects the axial force and the resulting equations, in terms of the displacements, are coupled nonlinear with variable coefficients. The analysis is performed with and without considering the axial deformation. This requires the solution of two different initial boundary value problems. The deviations of the two approaches are studied and compared with the linear theory. The solution of the problems was achieved using the analog equation method (AEM) of Katsikadelis [START_REF] Katsikadelis | The analog equation method. A boundary-only integral equation method for nonlinear static and dynamic problems in general bodies[END_REF] as it was developed for the nonlinear dynamic analysis of beams [START_REF] Katsikadelis | Non-linear dynamic analysis of beams with variable stiffness[END_REF]. According to this method, the two coupled nonlinear hyperbolic partial differential equations with variable coefficients are replaced by two uncoupled linear ones pertaining to the axial and transverse deformation of a substitute beam with unit axial and bending stiffness, respectively, under fictitious time dependent load distributions. Besides its accuracy, this method overcomes the shortcoming of a FEM solution, which experiences lack of convergence [START_REF] Andersen | Post-critical behaviour of Beck's column with a tip mass[END_REF] and depends on discretization. Example problems of uniform Beck' column and Beck' column with linearly varying height are presented, which illustrate the effectiveness of the employed method to handle this problem.

Moreover, useful conclusions are drawn concerning the influence of the axial deformation.

Thus, in case of uniform Beck' column this influence is negligible, while in case of Beck's column with linearly varying height is dominant and should be always included in the analysis.

A c c e p t e d

m a n u s c r i p t 4

Governing equations

Consider an initially straight beam of length l of viscoelastic material vibrating in a viscous medium. The beam has variable cross-section = ( ) A A x and moment of inertia ( ) = I I x . The x axis coincides with the neutral axis of the beam, which is bent in its plane of symmetry xz under the action of a tangential compressive follower tip load P (see Fig. 1a). We assume that there is no abrupt variation in cross-section of the beam so that the Euler-Bernoulli theory is valid [START_REF] Katsikadelis | Buckling load optimization of beams[END_REF]. In the following the equations are derived (a) for nonlinear analysis including the axial deformation, (b) for nonlinear analysis excluding the axial deformation and (c) for linear analysis.

Nonlinear theory including axial deformation

Moderate large deflections are considered. In this case the nonlinear kinematic relation retains the square of the slope of the deflection, while the strain component remains still small compared with the unity. Thus we have

ε κ = + + 2 1 2 ( , )
, ,

x x x x z u w z ( 1
)
where = ( , ) u u x t and = ( , ) w w x t are displacements along the x and z axis, respectively, and κ is the curvature of the deflected axis given as ( )

3/2 2 , 1 , xx x w w κ = - + (2)
The equations of motion are derived by considering the equilibrium of the deformed element. Thus, referring to Fig. 1b and taking into account the inertia and external damping forces we obtain

( ) cos sin , 0 - -+ - = x mu cu N Q θ θ (3) ( ) sin cos , 0 - - + + = x mw cw N Q θ θ (4)
, xx w κ = - (8)
The stress resultants are evaluated by integrating appropriately the normal stress

E E E E t t ε σ ε ε * *   ∂ ∂   = + = +       ∂ ∂ (9)
Thus, the axial force and the bending moment are obtained as [START_REF] Bolotin | Effect of damping on the postcritical behaviour of autonomous non-conservative systems[END_REF] where E * is the coefficient of dynamic visco-elastic resistance [START_REF] Hudson | The Excitation and Propagation of Elastic Waves[END_REF].

( ) ( ) ( ) 2 1 2 2 1 2 , , , , , 2 , , * *
  ∂   = + +       ∂ = + + + x x x x x x x N EA E A u w t EA u w E A u w w (10) , , , * *   ∂   = - +       ∂ = - - xx xx xx M EI E I w t EIw E Iw
Substituting Eqs. [START_REF] Kounadis | Hamiltonian weakly damped autonomous systems exhibiting periodic attractors[END_REF] into Eqs. ( 3) and (4) and using Eq. ( 5) to eliminate Q , we obtain equations of motion in the form

( ) , , , , 0 - -+ - = x x x x mu cu N M w (12) ( ) , , , 0 - - + + = xx x x mw cw M Nw ( 13
)
which by virtue of Eqs. [START_REF] Kounadis | Non-potential dissipative systems exhibiting periodic attractors in region of divergence[END_REF] and [START_REF] Bolotin | Effect of damping on the postcritical behaviour of autonomous non-conservative systems[END_REF] become

( ) ( ) ( ) 2 1 2 , , , 2 , , , , , , , , 0 * *
  - -+ + + +       + + =   x x x x x x xx xx x x x mu cu EA u w E A u w w EIw E Iw w (14)
A c c e p t e d m a n u s c r i p t

( ) ( ) ( ) { } 2 1 2 , , , , , , 2 , , , , 0 * *
- - - + +   + + + =     xx xx xx x x x x x x x mw cw EIw E Iw EA u w E A u w w w (15)
The pertinent boundary conditions are

( ) 0 0 u = and ( ) 2 1 2 , , x x x l N l P EA E A u w t * =   ∂     = -= + +         ∂ (16a,b) ( ) 0 0 w = and ( ) 0 , , xx x x l Q l EI E I w t * =     ∂     = = - +         ∂   (17a,b) ( ) , 0 0 x w = and ( ) 0 , xx x l M l EI E I w t * =   ∂   = = - +       ∂ (18a,b)
and the initial conditions are

( ) , 0 ( ) u x u x = , ( ) , 0 ( ) u x u x = (19a,b) ( ) , 0 ( ) w x w x = , ( ) , 0 ( ) w x w x = (20a,b)
where ( ), ( ), ( ), ( ) u x u x w x w x are prescribed spatial functions.

Without restricting the generality in our analysis we neglect the axial inertia and damping forces, whose influence will be the subject of further investigation. Thus the equations of motion are reduced to

( ) ( ) 2 1 2 , , , , , , , 0 
   + + =       x x x x x x x x EA u w EIw w (21) ( ) ( ) { } 2 1 2 , , , , , , , 0 * 
 - -- + + + =     xx xx xx x x x x mw cw EIw E Iw EA u w w (22)
Moreover, the boundary conditions after dropping the time dependent terms [START_REF] Kar | Stability of a nonuniform cantilever subjected to dissipative and nonconservative force[END_REF][START_REF] Rao | Stability of tapered cantilever columns subjected to a tipconcentrated follower force with or without damping[END_REF][START_REF] Claudon | Détermination et maximisation de la charge critique d' une colonne de Hauger en présence d' amortisement[END_REF] become

( ) 0 0 u = and
( )

2 1 2 , , x x x l EA u w P = + =- (23a,b) ( ) 0 0 w = and ( ) , , 0 xx x x l EIw = - = (24a,b) ( ) , 0 0 x w = and , 0 xx x l EIw = - = (25a,b)
while the initial conditions are limited only to Eqs. [START_REF] Hudson | The Excitation and Propagation of Elastic Waves[END_REF].

A c c e p t e d m a n u s c r i p t 7 2.2 Nonlinear theory excluding axial deformation

We start with Eq. ( 12), in which we drop the axial inertia and damping forces. This results in equation

( ) , , , , 0 - = x x x x N M w ( 26
)
which can be readily integrated independently to yield

, , = + x x N M w c ( 27
)
Then using the boundary conditions (16b) and (18b) for the beam end x l = , we obtain

= - c P (28)
Thus, the axial force is given as

( ) , , = --x x N x P M w (29)
Introducing the bending moment M from Eq. ( 11) without the contribution of the time dependent term into Eq. (29), we obtain ( )

( ) , , , = -- xx x x N x P EIw w ( 30
)
which is substituted into Eq. ( 13) to yield the counterpart of Eq. ( 15), when the axial deformation is neglected. Thus we obtain the nonlinear equation of motion excluding the axial deformation as

( ) ( ) 2 , , , , , , , , 0 xx xx xx xx xx x x x mw cw EIw E Iw Pw EIw w *   + + + + + =   (31)
The boundary conditions ( 24) and ( 25) hold also in this case together with the initial conditions (20).

Linear theory

In this case Eqs. ( 30) and (31) are simplified to [START_REF] Beck | Knicklast des einseitig eingespannten, tangential gedrücten Stabes[END_REF] A c c e p t e d m a n u s c r i p t

( ) = - N x P (32) ( ) , , , , 0 * + + + + = xx xx xx xx mw cw EIw E Iw Pw (33)
under the boundary conditions ( 24) and ( 25) and the initial conditions [START_REF] Hudson | The Excitation and Propagation of Elastic Waves[END_REF].

The AEM solution for the nonlinear dynamic analysis of Beck's column

Eqs. ( 21) and ( 22) are solved using the AEM, which for the problem at hand is applied as follows. Let ( , ) u u x t = and ( , ) w w x t = be the sought solutions, which are two and four times differentiable in ( ) 0,l , respectively. Noting that Eqs. ( 21) and ( 22) are of the second order with respect to u and of fourth order with respect to w , respectively, we obtain by differentiating 21) and [START_REF] Beck | Knicklast des einseitig eingespannten, tangential gedrücten Stabes[END_REF]. The fictitious loads are established by developing a procedure based on the integral equation method for one-dimensional problems.

( ) 1 , , xx u b x t = (34) ( ) 2 ,
Thus, the integral representations of the solutions of Eqs. (34) and (35) are written as ( ) ( )

1 2 1 1 0 ( , ) , , l u x t c x c G x b t d ξ ξ ξ = + + ∫ (36) ( ) ( )
G x ξ = - (38)
( )

2 1 12 2 G x x ξ ξ = - - (39)
are the fundamental solutions (free space Green's functions) of Eqs. ( 34) and (35), respectively.

The derivatives of u and w are obtained by direct differentiation of Eqs. ( 36) and (37).

This yields ( ) ( ) ( ) The interval ( ) 0,l is divided into N equal elements (see Fig. 2) on which 1 b and 2 b are assumed to vary according to a certain law (constant, linear, parabolic etc). The constant element assumption is employed here, because the numerical implementation becomes very simple and the obtained numerical results are very good.

1 1 1 0 , , , , , l x x u x t c G x b t d ξ ξ ξ = + ∫ , ( ) ( ) 1 , , , xx u x t b x t = (40a,b) ( ) ( ) ( ) 2 3 4 5 2 2 0 , , 3 2 , , , l x
x w x t c x c x c G x b t d ξ ξ ξ = + + + ∫ (41a) ( ) ( ) ( ) 3 4 2 2 0 , , 6 2 , , , l xx xx
w x t c x c G x b t d ξ ξ ξ = + + ∫ (41b) ( ) ( ) ( ) 3 2 2 0 , , 6 , , ,
After discretization of Eqs. (36) and (37) we obtain ( , ) ( ) ( ) ( , )

2 2 1 1 1 1 ( , ) ( ) ( ) ( , ) N j j j j j j u x t x c t b G x d ξ ξ ξ - = = = + ∑ ∑ ∫ (42)
N j j j j j j w x t x c t b G x d ξ ξ ξ - + = = = + ∑ ∑ ∫ (43) or 1 1 1 1 ( , ) ( ) ( ) u x t x x = + H c G b ( 44
) 2 2 2 2 ( , ) ( ) ( ) w x t x x = + H c G b (45)
where 1 ( ) x G and 2 ( ) x G are 1 N × known matrices originating from the integration of the kernels 1 ( , ) G x ξ and 2 ( , ) G x ξ on the elements, respectively;

1 ( )

1 x x   =     H and 3 2 2 ( ) 1 x x x x   =     H ; 1 1 2 { , } T c c = c ; 2 3 4 5 6 { , , , } T c c c c = c
; 1 b , 2 b are the vectors containing the values of the fictitious loads at the nodal points, respectively. Similarly, we obtain for Eqs. ( 40) and ( 41)

1 1 1 1 , ( , ) ( ) ( )
x x x u x t x x = + H c G b , 1 , ( , ) xx u x t = b (46a,b) 2 2 2 2
, ( , ) ( ) ( ) Finally, collocating Eqs. ( 21) and ( 22) at the N nodal points and substituting the relevant derivatives from Eqs. ( 46) and (47) yields the following equations of motion ( , , , )

x x x w x t x x = + H c G b , 2 2 2 2 , ( , ) ( ) ( )
xx xx xx w x t x x = + H c G b (47a,b)
1 1 2 1 2 ( , , , ) = K b b c c 0, (48)
+ - = Mb Cb K b b c c 0 , (49)

A c c e p t e d m a n u s c r i p t

where M and C are known N N × generalized mass and damping matrices and ) generalized stiffness vectors. The boundary conditions ( 23)- [START_REF] Katsikadelis | Boundary Elements: Theory and Applications[END_REF] are in general nonlinear and can be written as

1 2 1 2 ( , , , ) = f b b c c 0 (50)
Eq. ( 49) is the semi-discretized equation of motion of the Beck's column. The associated initial conditions result from Eq. ( 45) when applied to all nodal points and combined with Eqs. [START_REF] Hudson | The Excitation and Propagation of Elastic Waves[END_REF]. Thus, we have

( ) ( ) 1 2 2 22 0 - = - b G w Hc (51) ( ) 1 2 2 0 - = b G w (52)
The time step integration method for nonlinear equations of motion can be employed to solve Eq. ( 49). In each iteration for 2 b within a time step, the current value of 2 b is utilized to update the vectors 1 b and 1 2 , c c on the basis of Eqs. ( 48) and (50). This demands the solution of a nonlinear system of algebraic equations, which is performed using the modified Newton-Raphson method. In this paper, the average acceleration time step integration method was employed to solve Eq. ( 49) and the results were cross-checked by a time step integration method based on the analog equation method [START_REF] Katsikadelis | A new time step integration scheme for structural dynamics based on the analog equation method[END_REF]. Once the vectors

Numerical examples

On the base of the procedure described in previous section a FORTRAN program has been written for establishing the dynamic nonlinear response of damped, both externally and internally, beam columns with variable mass and stiffness properties subjected to follower . The employed initial conditions are ( )

4 3 2 ( ) 8 4 6 /24 l w x w ξ ξ ξ = - + and
( ) 0 w x = with / x l ξ =
and l w is the initial tip deflection of the beam. In Table 1 the computed critical loads are presented from (i) linear theory, (ii) nonlinear theory excluding the axial deformation and (iii) nonlinear theory including the axial deformation. The convergence of the method is shown by increasing the number of elements N . In all three theories, the computed values of the critical load are identical and coincide with those obtained by Andersen and Thomsen [START_REF] Andersen | Post-critical behaviour of Beck's column with a tip mass[END_REF], who used a perturbation analysis of the nonlinear equations of motion ignoring, however, the axial deformation of the beam. The FEM solution using geometric nonlinear dynamic finite element model for case (ii) gives results which are qualitative and quantitative the same, which, however, exhibit lack of convergence in some cases [START_REF] Andersen | Post-critical behaviour of Beck's column with a tip mass[END_REF].

In Table 2 the finite tip amplitude is presented from the nonlinear theory (i) excluding and (ii) including the axial deformation for various values of the axial load P . The two nonlinear theories give the same results. However, as compared with those obtained in [START_REF] Andersen | Post-critical behaviour of Beck's column with a tip mass[END_REF] only the case for 13.53 = P are found in good agreement. From the obtained results, it becomes apparent that the influence of the axial deformation on the critical load is negligible for uniform Beck's column. This is an interesting finding.

A c c e p t e d m a n u s c r i p t

13 Moreover, the time histories of the tip deflection for two different initial tip deflections l w are shown in Fig. 3 validating the supercritical bifurcation character of uniform Beck's column, while in Fig. 4 and Fig. 5 the variation of the critical load in regard to the external and internal damping is depicted, respectively. From Fig. 4, it can be pointed out that the critical load increases monotonically with the coefficient of the external viscous damping. On the contrary, the curve of Fig. 5 shows that the value of the critical load decreases to a minimum value ( 12.03

cr P = at 0.03 E E * =
) and thereafter it increases with increasing value of the coefficient of dynamic visco-elastic resistance (internal damping).

Example 2: Beck's column with variable mass and stiffness properties

The nonlinear dynamic stability of Beck's column with variable cross section has been studied. The employed data are the same with those in previous example. The height of the beam varies according to the linear law ( )

( ) 0 / 2 h x h a x l = + - (see Fig. 6) with 2 tan a φ =
being the taper ratio and 0 h the height at the half length. In order to compare the results with those of the previous example, the volume of the beam i.e.

0

V bh l = , was kept constant. The resulting beam should have no abrupt change of the cross-section so that the Euler-Bernoulli theory remains valid [START_REF] Katsikadelis | Buckling load optimization of beams[END_REF]. Boley [START_REF] Boley | On the accuracy of the Bernoulli-Euler theory for beams of variable section[END_REF] has shown that, for a beam with unit constant width, a rate of change of the cross-section 0.35 a yields an error of 7.5% , while for 0.17 a the error is 1.8% . This was also verified by the authors who treated the beam as a 2D elasticity problem and used the BEM to obtain the solution [START_REF] Katsikadelis | Boundary Elements: Theory and Applications[END_REF].

The computed critical loads from (i) linear theory, (ii) nonlinear theory excluding and (iii) including the axial deformation for various values of the ratio a are shown in Table 3.

The results obtained on the basis of the first two theories are coincident but differ considerably from those of the third one. This significant finding demands that the axial deformation on the nonlinear dynamic stability of Beck's column should be always included

Conclusions

In this paper the nonlinear dynamic stability of Beck's column with variable mass and stiffness properties in the presence of damping (both internal and external) has been investigated using a complete nonlinear dynamic analysis which includes the axial deformation. The solution of the derived coupled nonlinear equations of motion was achieved effectively using the analog equation method. This investigation has reached to certain striking effects concerning the influence of the axial deformation on the critical load of Beck's column with variable cross section. The main conclusions can be summarized as

A c c e p t e d m a n u s c r i p t

• The bifurcation is always supercritical in uniform and nonuniform Beck's column.

• The axial deformation affects considerably the critical load, when the beam has variable cross section. It may give remarkably lower or larger critical loads. Therefore, it should be always included in the analysis. On the contrary, it has negligible influence for uniform cross section.

• The exclusion of the axial deformation in the nonlinear dynamic analysis of Beck's column with variable mass and stiffness properties may yield larger loads, reducing thus the safety of the structure.

• The critical load increases monotonically in uniform and nonuniform Beck's column with the coefficient of the external viscous damping.

• The critical load decreases to a minimum value and thereafter increases for further increase of the value of the coefficient of dynamic visco-elastic resistance (internal damping) in uniform and nonuniform Beck's column.

• In the nonlinear theory including the axial deformation, the critical load increases monotonically with the taper ratio. This suggests shifting of the material towards the free end in order to maximize the critical load

• The shape of the bifurcation branches becomes steeper as the taper ratio decreases.

  density per unit length and c is the coefficient of the external viscous damping. For the case of moderate large deflections the following relations are valid[START_REF] Katsikadelis | Non-linear dynamic analysis of beams with variable stiffness[END_REF]

  34) and (35) describe the axial and bending linear response of a beam with constant unit axial and flexural stiffness subjected to the fictitious time dependent axial 1 b , and transverse 2 b , respectively. They indicate that the solution of Eqs. (21) and (22) can be established by solving Eqs. (34) and (35) under the boundary conditions (23)-(25), provided that the fictitious load distributions 1 b , 2 b are first determined. Eqs. (34) and (35) are quasistatic, that is the time is considered as a parameter. Note that Eqs. (34) and (35) are referred to as the analog equations to Eqs. (

  derivatives into Eqs. (21) and (22) yields the equations, from which the fictitious sources 1 b and 2 b can be determined. This can be implemented only numerically as follows.

A c c e p t e d m a n u s c r i p t 10

 10

  t =and their derivatives at any instant t are evaluated from Eqs. (44) through (47).

A c c e p t e d m a n u s c r i p t

  forces. The uniform and nonuniform Beck's column has been studied as an illustrative example of the developed analysis and solution method.Example 1: Uniform Beck's columnThe nonlinear dynamic stability of Beck's column with uniform rectangular cross-

A c c e p t e d m a n u s c r i p t 14 in

 14 the analysis when the beam has variable mass and stiffness properties. It should be also observed in the nonlinear theory (iii) that the critical load increases monotonically with the taper ratio, which means that the material must be shifted towards the free end in order to obtain the maximum critical load.Moreover, in Fig.7are shown the bifurcation diagrams from nonlinear theory (iii) for various values of the taper ratio. The shape of the bifurcation branches becomes steeper as the taper ratio decreases. For this reason the case of taper ratio 0.15 a =is investigated for possible subcritical bifurcation. Fig.8and Fig.9show the variation of the critical load in regard to the external and internal damping, respectively. It is apparent from these figures that the same qualitatively conclusions with those of uniform Beck's column can be drawn.Namely, the critical load increases monotonically with the coefficient of the external viscous damping, while it is decreasing to a minimum value ( further increase of the internal damping. Finally, in order to establish the bifurcation character, two nonlinear dynamic analysis are performed, one below ( load. The time histories of the tip deflection are shown in Fig. 10 validating the supercritical bifurcation character of nonuniform Beck' column.

A c c e p t e d m a n u s c r i p t A c c e p t e d m a n u s c r i p t

Fig. 1 .Fig. 2 .

 12 Fig. 1. (a) Beck' column with variable mass and stiffness properties and (b) forces and moments acting on the deformed element.

Fig. 3 .

 3 Fig. 3. Example 1: Time history of the tip deflection of uniform Beck' column ( 0.1 c = , 0.01 E E * = , 12.94 cr P = ) for 13.33 P = (i) 0.05 l w = and (ii) 0.20 l w = .

A c c e p t e d m a n u s c r i p t 20 Fig. 4 .

 204 Fig. 4. Example 1: Critical load of uniform Beck' column versus external damping ( 0.01 E E * = , 30 = N).

Fig. 5 .

 5 Fig. 5. Example 1: Critical load of uniform Beck' column versus internal damping ( 0.1 = c , 30 = N).

Fig. 6 .

 6 Fig. 6. Linear variation of the height of the Beck' column in Example 2.

A c c e p t e d m a n u s c r i p t 21 Fig. 7 .

 217 Fig. 7. Example 2: Bifurcation diagrams of Beck' column with linearly varying height from nonlinear theory including the axial deformation for various values of the taper ratio a ( 0.1 c = , 0.01 E E * = , 30 = N).

Fig. 8 .

 8 Fig. 8. Example 2: Critical load of Beck' column with linearly varying height ( 0.15 a = -) versus external damping ( 0.05 E E * = , 30 = N).

A c c e p t e d m a n u s c r i p t 22 Fig. 9 .Fig. 10 .

 22910 Fig. 9. Example 2: Critical load of Beck' column with linearly varying height ( 0.15 a = -) versus internal damping ( 0.1 = c , 30 = N).

Table 1 .

 1 Example 1: Critical load of Beck's column (

	c =	0.1	,	E	* =	0.01 E	) from (i) linear

Table 2 .

 2 Example 1: Finite tip amplitude of Beck's column (

	c =	0.1	,	E	* =	0.01 E	) from

Table 3 .

 3 Example 2: Critical load of Beck's column with linearly varying height (

	c =	0.1	,

A c c e p t e d m a n u s c r i p t


Angelo Di Egidio

Angelo Luongo

Achille Paolone

Linear and nonlinear interactions between static and dynamic bifurcations of damped planar beams

Keywords: Stability analysis, damping effects on stability, non-conservative loads, bifurcation, Multiple Scales Method, direct perturbation approach, divergence and Hopf bifurcations, beam continuous model

published or not. The documents may come

Linear and nonlinear interactions between static and dynamic bifurcations of damped planar beams

Columns subjected to follower forces, after the pioneering paper by Beck [START_REF] Beck | Die Knicklast des einseiting eingespannten tangential gedrückten Stabes[END_REF], have recently attracted the attention of many researchers, particularly in aerospace, where tangential forces are produced by jets and rocket motors [2]. Columns with a tip mass have also been studied, to account for the influence of the mass of the device delivering the load [START_REF] Andersen | Post-critical behavior of Beck's column with a tip mass[END_REF]. Lumped damping forces and couples have then been considered to model the external action of air on the mass [START_REF] Detinko | Lumped damping and stability of Beck column with a tip mass[END_REF]. On the other hand, passive control systems for vibration reduction also demand the introduction of dissipation devices (dashpots), which can remarkable affect the critical and postcritical column behavior. The effect of dashpots on the linear stability of the Beck column has, for example, been studied in [START_REF] Kirillov | Dissipation induced instabilities in continuous non-conservative systems[END_REF]; a general treatment of the effect of an added small damping on the stability of linear continuous nonconservative systems can be found in [START_REF] Kirillov | The effect of small internal and external damping on the stability of distributed non-conservative systems[END_REF]. However, except for [START_REF] Andersen | Post-critical behavior of Beck's column with a tip mass[END_REF], where the Multiple Scale Method has been employed to investigate supercritical and subcritical Hopf bifurcations, the attention of the aforementioned researchers has been limited to the linear stability problem. Here, it is believed important to extend the analysis to the nonlinear range, accounting for possible interactions between static and dynamic bifurcation mechanism. Such a problem has been deeply studied in [START_REF] Kounadis | On the Failure of Static Stability Analyses of Nonconservative Systems in Regions of Divergence Instability[END_REF]- [START_REF] Kounadis | Non-Potential Dissipative Systems Exhibiting Periodic Attra-ctors in Region of Divergence[END_REF],

where some discrete versions of the continuous column have been considered.

It is well-known that the essential dynamics of a nonlinear system, finite-or infinitedimensional, close to a bifurcation point, is captured by an equivalent, lower dimensional, reduced system [START_REF] Arnold | Geometrical Methods in the Theory of Ordinary Differential Equations[END_REF]- [START_REF] Troger | Nonlinear Stability and Bifurcation Theory[END_REF]. Reduction method, as linear and nonlinear Galerkin, center manifold and approximate inertial manifold methods, have been thoroughly discussed in [START_REF] Steindl | Methods for Dimension Reduction and Their Application in Nonlinear Dynamics[END_REF]. In the last decade, as an alternative and more engineering oriented approach, the authors have systematically applied the Multiple Scales Method to analyze a number of bifurcations of linear codimension-one, two and three, to general finite dimensional systems [START_REF] Luongo | Perturbation Methods for Bifurcation Analysis from Multiple Nonresonant Complex Eigenvalues[END_REF]- [START_REF] Luongo | Multiscale analysis of defective multiple-Hopf bifurcations[END_REF]; a review paper [START_REF] Luongo | Multiple scale bifurcation analysis for finitedimensional autonomous systems[END_REF] resumes their main results. More recently, they have extended the method to infinite dimensional system, to analyze divergence, Hopf and double-zero bifurcations [START_REF] Luongo | Divergence, Hopf and double-zero bifurcations of a nonlinear planar beam[END_REF]- [START_REF] Paolone | Flexural-torsional bifurcations of a cantilever beam under potential and circulatory forces: II Post-critical analysis[END_REF]. The method is based on the direct treatment of the original (integro)-differential equations, avoiding any a priori discretization (as that, e.g., performed in [START_REF] Andersen | Post-critical behavior of Beck's column with a tip mass[END_REF]), according

A c c e p t e d m a n u s c r i p t

to the so-called direct method, widely applied by Nayfeh and co-workers [START_REF] Nayfeh | Reduced-Order Models of Weakly Nonlinear Spatially Continuous Systems[END_REF]- [START_REF] Rega | Multiple Resonances in Suspended Cables: Direct versus Reduced-Order Models[END_REF], and many other authors (see, e.g., [START_REF] Lacarbonara | Elastodynamics of nonshallow suspended cables: Nonlinear vibration characteristics[END_REF][START_REF] Lacarbonara | Galloping instabilities of geometrically nonlinear nonshallow cables under steady wind flows[END_REF]), to several problems of nonlinear dynamics.

In this paper the algorithm is applied to study a cantilever beam, constrained by a spring and two dashpot, loaded by a follower force. This system, for its simplicity, was already studied in [START_REF] Luongo | Divergence, Hopf and double-zero bifurcations of a nonlinear planar beam[END_REF][START_REF] Luongo | Bifurcation equations through Multiple-Scales analysis for a continuous model of a planar beam[END_REF] as an example of structure undergoing divergence, Hopf and double-zero bifurcation. A deeper parametric analysis performed on the system permitted to reveal a richer bifurcation scenario, including Hopf-divergence and resonant and non-resonant double-Hopf, not discovered in the previous analysis. Therefore, the beam viscous-elastically restrained could be taken as paradigmatic system undergoing all the low-codimension bifurcations of mechanical interest. Here, the analysis of the codimension-two bifurcations is completed, while the codimension-three one (resonant Hopf-Hopf) is left for future investigation.

The paper is thus organized. In Sect. 2 the equations of motions are given, and detailed in Appendix A, where they are derived through a procedure alternative to that of Ref. [START_REF] Luongo | Divergence, Hopf and double-zero bifurcations of a nonlinear planar beam[END_REF].

Moreover the linear adjoint problem is defined. In Sect. 3 the critical scenario is depicted.

Considerably emphasis is given to the influence of the relative magnitude of the two dashpots on bifurcations. The well-known destabilizing effect of damping is detected, and a new paradoxical result discovered. In Sect. 4 the post-critical analysis is carried-out. The bifurcation equations are first derived and then numerically studied to built-up equilibrium paths and bifurcation diagrams. Finally, in Sect. 5 some conclusions are drawn.

Model

A planar beam is considered, fixed at end A and constrained by a linear visco-elastic device at end B, loaded at the tip by a follower force of intensity P (Fig. 1). The device consists of an extensional spring of stiffness k e and two dashpots of constants c e and c t , of extensional and torsional type, respectively. The beam is assumed to be inextensible and shear-undeformable, so that bending is the unique strain measure.

A c c e p t e d m a n u s c r i p t 2.1 The equations of motion

The equations of motion of the beam were derived in [START_REF] Luongo | Divergence, Hopf and double-zero bifurcations of a nonlinear planar beam[END_REF] through a variational approach, by eliminating the rotation through the no-shear constraint and by using a lagrangian multiplier to account for the axial constraint. Then, after condensation of the longitudinal displacement and of the lagrangian multiplier, the following, unique, nonlinear integro-differential equation, corrected up-to the third order, was obtained in nondimensional form:

ü + u + 2µu + u 3 + 4u u u + u 2 u + µ 3u 2 -u 2 B u (1)
+ s 1 s 0 1 2 u 2 ds .. ds u + s 0 1 2 u 2 ds .. u = 0
with the relevant boundary conditions:

-u B + κu B + ξ e uB = -µu 3 B + u B u 2 B + u 2 B u B , u B + ξ t u B = -1 2 ξ t u 2 B u B -1 2 u B u 2 B u A = 0, u A = 0 (2)
where A and B denote evaluation at the beam ends. In Eqs. ( 1)-( 2) the following nondimensional variables and parameters were introduced:

t = ωt, s = s/l, ũ = u/l, µ = P l 2 /2EI κ = k e l 3 /EI, ξ e = c e ωl 3 /EI, ξ t = c t ωl/EI (3)
where u(s, t) ≡ u y (s, t) is the transversal displacement of the beam at the abscissa s and time t (with tilde dropped), l is the length, ω = (EI/ml 4 ) 1/2 is a frequency and EI the bending stiffness. Moreover, µ is the load parameter, κ is a nondimensional stiffness parameter and ξ e and ξ t two nondimensional damping coefficients, referred ahead as the 'spring stiffness' and the 'extensional and torsional damping', respectively.

To make the paper self-contained, the equations of motion (1) and (2) are reobtained in Appendix A by following a direct method, based on equilibrium and later condensation of all the internal reactive forces (axial and shear), which are not achievable by a variational procedure. The direct method, moreover, reveals the existence of an extra-term, lost in the previous formulation (namely the nonlinearity affected by µ). The slight difference in the two derivations was due to the (usually accepted) truncation at the second-order, instead of

A c c e p t e d m a n u s c r i p t

the third-order, of the axial strain-displacement relation used in the variational principle. A recomputation of the variational procedure, accounting for this higher-order term, permitted to recover the results of the direct method.

The linear direct and adjoint problems

The linearized problem Eqs. ( 1)-( 2) is first studied. By letting u(s) = φ(s) exp(λt), the following eigenvalue problem is obtained:

φ + 2µφ + λ 2 φ = 0 φ A = 0, φ A = 0 -φ B + κφ B + λξ e φ B = 0, φ B + λξ t φ B = 0 (4)
In Eq. ( 4), λ ∈ C is the eigenvalue and φ(s) ∈ C is the associated (right) eigenvector.

It is useful to introduce for later convenience also the adjoint eigenvalue problem. It is derived from Eq. ( 4) performing scalar multiplication (in the complex field) by the dual variables ψ(s) ∈ C and, successively, integrating by parts (i.e. enforcing the so-called bilinearidentity). The following equations and boundary conditions was obtained:

ψ + 2µψ + λ2 ψ = 0 ψ A = 0, ψ A = 0 -ψ B + κψ B -2µψ B + λξ e φ B = 0, ψ B + 2µψ B + λξ t ψ B = 0, (5)
where the overbar denotes complex conjugate. In Eqs. (5), ψ(s) is the left eigenvector associated with the eigenvalue λ (the latter being identical in the two problems).

It should be noted that the two field equations (4 1 ) and (5 1 ) are identical, while the mechanical boundary conditions (4 4,5 ) and (5 4,5 ) differ each other; this is a consequence of the fact the nonconservative sources (follower force and dashpots) are located at the end of the beam.

A c c e p t e d m a n u s c r i p t 3 Linear Stability Analysis

Stability of the trivial, rectilinear configuration, of the beam is governed by the eigenvalue problem Eq. ( 4). It admits the eigenvector:

φ(s) = c 1 (cos ps -cosh qs) + c 2 (sin ps - p q sinh qs) (6)
where:

q 2 := µ 2 -λ 2 -µ, p 2 := µ 2 -λ 2 + µ (7)
have been set. The arbitrary constants c = (c 1 , c 2 ) T , together with the eigenvalue λ, satisfy the algebraic equations:

S(λ) c = 0 (8)
where:

S(λ) =         -(κ + λ ξ e ) (cos p -cosh q)
+p 3 sin p -q 3 sinh q, -(κ + λ ξ e ) (sin p -p/q sinh q)

-p 3 cos p -q 2 p cosh q -λ ξ t (p sin p + q sinh q)

-p 2 cos p -q 2 cosh q, λ ξ t p(cos p -cosh q)

-p 2 sin p -p q sinh q         (9)
is the 'dynamical stiffness matrix' of the system, depending on the eigenvalue λ, in addition to the system parameter, µ, κ, ξ e and ξ t .

In order to find divergence and Hopf boundaries in the parameter space, the loci of the roots λ = 0 and λ = iω of the characteristic equation detS(λ) = 0 must be found, respectively.

In this analysis, the load µ and the stiffness κ are taken as bifurcation parameters, while the damping coefficients ξ e and ξ t are assumed to be auxiliary parameters. In other words, sections ξ e = const, ξ t = const, of the four-dimensional parameter space are considered.

Divergence occurs at the manifold D on which detS(0) = 0. Since this equation presents an undetermined form, a limit process for λ → 0 has to be performed. By letting p → √ 2µ, q → iλ √ 2µ in detS(λ), expanding for small λ and retaining only the leading terms, the following equation is obtained:

2.82µ 3/2 κ cos( 2µ) -5.66µ 5/2 -2µκ sin( 2µ) = 0 (10)

A c c e p t e d m a n u s c r i p t

which, implicitly, defines a multi-branch curve D in the (κ, µ)-plane.

Hopf bifurcation occurs at the manifold H on which detS(iω) = 0; in addition, [∂λ/∂n] = 0 must hold, n being the normal to the curve in the plane (transversality condition). This last condition permits to exclude degenerate Hopf cases, as that ones studied by Kounadis [START_REF] Kounadis | Hamiltonian weakly damped autonomous systems exhibiting periodic attractors[END_REF]. When λ = i ω, p and q are real, so that the unique imaginary terms in S are those due to damping. By separating real and imaginary terms in the characteristic equation, two real equations of the type f i (κ, µ, ω; ξ e , ξ t ) = 0 (i = 1, 2) are found. They, again, implicitly define a multi-branch curve H in (κ, µ)-plane, parameterized in the ω-parameter. Transversality is either checked numerically, by evaluating [λ] closely to H, or analytically, via an eigenvalue sensitivity analysis (see Appendix B). It is worth noting that, in general, Hopf bifurcations depend on both ξ e and ξ t parameters. However, when only one of the damping coefficients is different from zero, the complex characteristic equation splits in two real equations of the type f (κ, µω) = 0, ξg(κ, µ, ω) = 0 (with ξ = ξ e or ξ = ξ t ), i.e. the Hopf curves are independent of the (unique) damping.

To better illustrate the system behavior, and emphasize its dependence on damping, some sub-systems, obtained by removing one of both dashpots from the complete system, are preliminarily studied. All these systems, of course, exhibit the same divergence curve, while the Hopf curves change. for large κ's, of the dynamic instability is due to the fact that the system tends to the clumped-support beam, for which the load is of conservative type (see [START_REF] Leipholz | Stability theory[END_REF]).

The previous scenario abruptly changes, if an extensional damper is added (Fig. 3).

Irrespective of the magnitude of ξ e (for the reason above mentioned), the Hopf curve goes down for all values of κ, and tends to tangentially merge with the lower part of curve D,

A c c e p t e d m a n u s c r i p t

similarly to known results [START_REF] Kirillov | Dissipation induced instabilities in continuous non-conservative systems[END_REF][START_REF] Hermann | Instability modes of canitilevered bars induced by fluid flow through attached pipes[END_REF]. Therefore, due to extensional damping, dynamic instability occurs also for large κ's. However, for small κ's, ξ e has a destabilizing effect, according to the well-known Ziegler paradox.

The effect of the torsional damping, when it acts alone (ξ e = 0), is then studied (Fig. 4).

First, the occurrence of Hopf bifurcations from zero values of the load µ, as suggested by the numerical analysis, were analytically checked. In this special µ = 0 case, p = q = √ ω (Eq. 7), follow. The characteristic equation detS(iω) = 0 then reduces to the very simple form sin 2 p = 0, which admits the double roots p n = nπ (n = 1, 2, . . .). Therefore, at special values κ n (n = 1, 2, . . .) of κ, and µ = 0, two Hopf-curves cross each-other, with the same values of the frequency, namely ω n = n 2 π 2 . Thus, a infinite sequence of 1 : 1 resonant double Hopf bifurcations (RH) takes place on the κ-axis. It is worth noting that, for such sub-system, dynamic instability manifests itself both for µ > 0 (compressive) and µ < 0 (tensile) loads. The Hopf curves depicted in Fig. 4 were numerically evaluated, and turned out to be independent of ξ t . Since they cross the divergence curve and, moreover, cross each other, they reveal the existence of non-resonant Hopf-Hopf points (HH) and Hopf-divergence points (HD), marked in the figure. In addition, a critical point HHD of higher codimension appears, in which two Hopf curves and a divergence curve cross each other, two of them being tangent.

An other interesting aspect offered by the dynamic bifurcations of the purely torsionally damped system is the following. At a Hopf point, since the first of the two Eqs. (8) (i.e. the 'shear' boundary condition) has real coefficients (rememeber that p and q are real and ξ e = 0), then c 1 and c 2 are real, and the eigenvector ( 7) is (quite unusually) also real. This entails that, in the second of Eqs. (8) (i.e. in the 'moment ' boundary condition), the bending moment (real) and the damping couple (imaginary) must vanish separately. Therefore, at bifurcation, the tip of the beam does not rotate and the follower force, paradoxically behaves as a conservative force. However, this exceptional circumstance does not hold any more close to the curve H, since perturbations of the parameters active the nonconservative mechanism.

Indeed, the transversality conditions turned out to be satisfied.

A c c e p t e d m a n u s c r i p t

As a further step, a small extensional damping was added to the purely torsionally damped system. First, it was analytically checked that, due to ξ e , the resonant double Hopf points vanish, i.e. ξ e is an imperfection parameter for such bifurcation. A numerical analysis led to the plot of Fig. 5. This appear to be a slight perturbation of the plot of Fig. 4, if only the µ > 0 half-plane is considered. Hopf points was instead not found in the µ < 0 half-plane.

The perturbation, moreover, makes the HHD point to vanish.

Comparable values of the two damping coefficients was finally consider in Figs. 6 and7. They show (Figs. 6 a and 7 a ) that large extensional damping does not qualitatively modify the scenario of Fig. 5 (small ξ e ). Therefore, it can be concluded, that the system behavior is more strongly affected by torsional damping; extensional damping acts as an imperfection, destroying some bifurcation points, but slightly modifying the curves. large ξ e -damping, the close Hopf curve does not intersect the divergence curve. The analysis, when repeated for a fixed κ and different values of the ξ t -parameter, leads to the results of Fig. 7 b . This shows that, when ξ e vanishes, the critical load is independent of ξ t , consistently with results of Fig. 4.

Post-Critical Analysis

The Multiple Scale Method [START_REF] Nayfeh | Introduction to Perturbation Techniques[END_REF] is applied to analyze the system behavior at a Hopf-divergence and at a non-resonant double-Hopf points. The fundamental steps of the method are resumed

A c c e p t e d m a n u s c r i p t

here. A perturbation parameter ε is introduced as a measure of the distance of a generic point from a bifurcation point in the parameter space (κ, µ). Several independent timescales t k := ε k t are defined and the state variable u(s) is expanded in a Taylor series of ε.

By equating terms of the same power of ε, linear perturbation equations having the same operator are obtained, and then solved in sequence for the series coefficients. Except for the lower-order eigenvalue problem, higher-order equations are non-homogeneous; they admit a solution if and only if the known term belongs to the range of the singular operator, i.e. if it is orthogonal to the solutions of the adjoint homogeneous problem [START_REF] Luongo | Perturbation Methods for Bifurcation Analysis from Multiple Nonresonant Complex Eigenvalues[END_REF]. Solvability equations furnish the bifurcation equations, governing the asymptotic dynamics of the system, reduced to the center manifold. The procedure is sketched ahead.

Bifurcation Equations

The stiffness κ and the load µ are taken as bifurcation parameters, and their increments with respect to the critical values κ 0 and µ 0 , denoted by, and ordered as:

β := κ -κ 0 = O(ε 2 ), γ := µ -µ 0 = O(ε 2 ) (11)
By introducing the series expansions:

u = εu 1 + ε 3 u 3 + . . . , d dt = d 0 + ε 2 d 2 + . . . (12)
with d k = ∂/∂t k and t k = ε k t (k = 0, 2, . . .), the following perturbation equations and boundary conditions are derived:

Order ε:                d 2 0 u 1 + u 1 + 2µ 0 u 1 = 0 u 1A = 0, u 1A = 0 ξ e d 0 u 1B -u 1B + κ 0 u 1B = 0 ξ t d 0 u 1B + u 1B = 0 (13)

A c c e p t e d m a n u s c r i p t

Order ε 3 :

                           d 2 0 u 3 + u 3 + 2µ 0 u 3 = -2d 0 d 2 u 1 -2γu 1 -u 3 1 -4u 1 u 1 u 1 -u 2 1 u 1 -µ 0 3u 2 1 -u 2 1B u 1 - s 1 d 2 0 s 0 1 2 u 2 1 ds ds u 1 - s 0 d 2 0 1 2 u 2 1 ds u 1 u 3A = 0, u 3A = 0 ξ e d 0 u 3B -u 3B + κ 0 u 3B = -ξ e d 2 u 1B -βu 1B -µ 0 u 3 1B + u 1B u 2 1B + u 2 1B u 1B ξ t d 0 u 3B + u 3B = -ξ t d 2 u 1B -1 2 ξ t u 2 1B d 0 u 1B -1 2 u 1B u 2 1B (14)
Equations ( 14) (generating solution) admit the following non-decaying solution:

u 1 = A 1 (t 2 , t 4 , . . .)φ 1 (s)e iω 1 t 0 + A 2 (t 2 , t 4 , . . .)φ 2 (s)e iω 2 t 0 + c.c. (15)
or

u 1 = a 1 (t 2 , t 4 , . . .)φ 1 (s) + A 2 (t 2 , t 4 , . . .)φ 2 (s)e iωt 0 + c.c. ( 16
)
for non-resonant double Hopf or Hopf-divergence bifurcation, respectively. In Eqs. (15) A j = 1/2 a j exp(iθ j ) (j = 1, 2) are complex constants with real amplitude a j and phase θ j , φ j (s) are complex eigenfunctions and c.c. stands for the complex conjugate of preceding terms.

The same notation holds for Eq. ( 16), with a 1 ∈ R and φ 1 (s) ∈ R. Substituting Eqs. ( 15) or [START_REF] Luongo | Multiple Time Scale Analysis for Bifurcation from a Multiple-Zero Eigenvalue[END_REF] in Eqs. (14 1 ) and enforcing solvability, bifurcation equations, formally identical for the two cases, are drawn:

   ȧ1 = (α 1β β + α 1γ γ) a 1 + R 122 a 1 a 2 2 + R 111 a 3 1 + O |a 1 | 5 + |a 2 | 5 ȧ2 = (α 2β β + α 2γ γ) a 2 + R 112 a 2 1 a 2 + R 222 a 3 2 + O |a 1 | 5 + |a 2 | 5 (17)
The coefficients α's and R's appearing in Eqs. [START_REF] Luongo | Multiscale analysis of defective multiple-Hopf bifurcations[END_REF] was obtained via a symbolic manipulation program [START_REF][END_REF]. Since they assume quite cumbersome forms, they are not reported here. However, the numerical values that they assume at the bifurcation points studied ahead (Sect.

4.2), are given in the Appendix C. Equations ( 17) constitute the bifurcation equations, in standard normal form and specialized to a symmetric system, for two pairs of nonresonant purely imaginary critical eigenvalues, or for a one zero and a purely imaginary pair of critical eigenvalues [START_REF] Troger | Nonlinear Stability and Bifurcation Theory[END_REF]. They are invariant under the transformations a 1 → -a 1 and/or a 2 → -a 2 .

A c c e p t e d m a n u s c r i p t

In order to build-up the bifurcation diagram, it is first necessary to evaluate the steadystate solutions of the dynamical system (17) and, then, to perform a stability analysis. The fixed points of Eqs. ( 17) are determined by setting ȧ1 = ȧ2 = 0. These solutions, in terms of the original system, correspond to two-frequency quasi-periodic motions around the rectilinear configuration of the beam or to one-frequency periodic motion around a nontrivial (buckled) position, respectively. Equations ( 17) admit the trivial solution a 1T = a 2T =0. Nontrivial steady-state solutions, with one or two non-vanishing components are then sought. If a 2 =0, equation Eq. ( 172 ) is identically satisfied, while Eq. ( 171 ) yields:

a 2 1 = - α 1β β + α 1γ γ R 111 (18)
Similarly, if a 1 =0, Eq. ( 171 ) is identically satisfied, while Eq. ( 172 ) yields:

a 2 2 = - α 2β β + α 2γ γ R 222 (19)
If both a 1 and a 2 are different from zero, Eqs. [START_REF] Luongo | Multiscale analysis of defective multiple-Hopf bifurcations[END_REF] give:

a 2 1 = (α 2β R 221 -α 1β R 222 )β + (α 2γ R 221 -α 1γ R 222 )γ R 111 R 222 -R 221 R 211 a 2 2 = (α 2β R 111 -α 1β R 211 )β + (α 2γ R 111 -α 1γ R 211 )γ R 221 R 211 -R 111 R 222 (20)
Since a 1 and a 2 are real, solutions ( 18)-( 20) exist only in a sub-domain of the control parameter plane. Finally, stability of the fixed points is analyzed by evaluating the eigenvalues of the Jacobian matrix at equilibria.

Numerical Results

The reduced dynamical system (17) is numerically analyzed for two sample systems, respectively at a Hopf-Hopf and at Hopf-divergence points. Figure 8 illustrates the HH and HD bifurcation points, the relevant critical parameters κ 0 and µ 0 .

The HH bifurcation is studied first. The bifurcation diagram is shown in Fig. 9. The parameter plane is divided in six sectors, delimited by four straight lines passing for the origin. Lines r 1 and r 2 , individually delimit the half-plane in which monomodal periodic

A c c e p t e d m a n u s c r i p t

solutions (a 1 , 0) and (0, a 2 ) respectively exist, according to Eq. ( 18) and ( 19); lines r 11 and r 22 , together, delimit the sector in which the bimodal quasi-periodic solution (a 1 , a 2 ) exists, according to Eqs. [START_REF] Luongo | Bifurcation equations through Multiple-Scales analysis for a continuous model of a planar beam[END_REF]. A sketch of the phase portraits relevant to the six sectors is also shown in Figure . It appears that stable solutions only exist in sectors 2 (trivial solution), 3, where the motion is periodic, and in sector 4, where the motion is bi-periodic.

Some equilibrium paths, denoted by I to IV in Fig. 9, are shown in Fig. 10. They can be thought as planar sections of three-dimensional plots (a i , β, γ), with i = 1, 2. Points A to L marked in Fig. 10 correspond to homonymous points in Fig. 9. It is seen that (Fig. 10 At point G un unstable a 1 -motion arises; at point H a stable a 2 -motion takes place. This exchange stability at I with the stable (a 1 , a 2 )-motion. Path IV is finally shown in Fig. 10 d , that displays as the a 1 -amplitude of the bi-periodic motion does not depend on the load

γ = µ -µ 0 .
The HD bifurcation is then addressed. The relevant bifurcation diagram is plotted in Fig. 11. Here, lines r 1 , r 2 , r 11 and r 22 keep the meaning previously illustrated, where the index 1 refers to divergence and index 2 to Hopf. As a peculiarity of this system, lines r 2 and r 22 are nearly undistinguishable; as a result, the parameter plane is divided in only five sectors. The relevant phase-portraits, in contrast with the previous case (Fig. 9), reveal the existence of a stable solution in the whole plane. The two equilibrium paths shown in Fig. 12 are therefore sufficient to illustrate the post-critical behavior. Along path I (Fig. 12

Conclusions

The critical and post-critical behavior of a beam constrained by lumped spring and dashpots, loaded by a follower force, has been studied. The following conclusions are drawn.

1. The beam undergoes Hopf, divergence, double-zero, Hopf-divergence and double-Hopf bifurcations, resonant or not, for different combinations of parameters.

2. If both dashpots are removed, an Hamiltonian Hopf bifurcation takes place along a curve tangent to the divergence curve, similarly to known results.

3. If only a dashpot of extensional type is added to the undamped beam, however small is the damping coefficient, it has a destabilization effect for any values of the spring stiffness, and the well-known Ziegler paradox takes place. Discontinuity exists between the undamped and the damped system.

4. If only a torsional damper is added to the undamped beam, an infinite number of resonant double-Hopf manifests itself for zero values of the load, so that simple Hopf bifurcations occur both for positive (compressive) and negative (tensile) forces. However, if a small extensional damping is added, such bifurcations are destroyed and instability for traction disappears.

5. In the purely torsional damping case, the critical mode at a generic Hopf point entails vanishing of both bending and rotation at the tip of the beam. Therefore, as a paradox, the nonconservative force behaves as a conservative one at the bifurcation.

6. A non degenerate double-zero bifurcation, occurring at the crossing of a Hopf and a divergence curve, where the first one dies, only can exist if both dampers are present in the structure.

A c c e p t e d m a n u s c r i p t

7. The scenario of the complete model can usefully be interpreted as a perturbation, due to the extensional damping, of the scenario of the purely torsionally damped structure.

8. The bifurcation equations relevant to the nonresonant double-Hopf and Hopf-divergence bifurcations are formally similar. They lead to a bifurcation scenario in a two-dimensional parameter space in which periodic motions around the trivial or buckled configurations, as well as bi-periodic motions, are present. Some sections of the amplitude-parameter surfaces have been given.

A Equations of Motion

A Lagrangian description of the motion is adopted (Fig. 1 b ). The beam is straight in its rest reference configuration C 0 , coincident with the prestressed equilibrium state, with directors D j , (j = x, y), oriented as the orthonormal versors of a fixed inertial reference frame (A; a j ).

The reference configuration is therefore described by X(s) = s a x where s ∈ (0, l) denotes an abscissa and l is the beam length. The referential description of the current beam configuration C t at time t is specified by the position vector x(s, t) = X(s) + u(s, t) with u(s, t)

the vector displacement, and by the proper orthogonal tensor R(s, t) describing the rotation that leads the directors D j to match the current directors d j (i.e. d j = RD j ); in the fixed inertial reference (A; a j ), it admits the following scalar representation:

R =   cos ϑ -sin ϑ sin ϑ cos ϑ   ( 21
)
where ϑ is the amplitude rotation around the axes oriented by the versor a z = a x × a y .

Denoting by (u x , u y ) the Cartesian components of u in the basis (a x , a y ), the current configuration C t is described by the longitudinal displacement u x (s, t), the transversal displacement field u y (s, t) and the rotation ϑ(s, t). The three displacements, however, are not independent, because of the internal constraints: 22), the remaining passive variables (u x , ϑ), truncated at the third order, are obtained:

tan ϑ = u y 1 + u x , (1 + u x ) 2 + u y 2 -1 = 0 (22)
u x = - 1 2 s 0 u 2 dx, ϑ = u + 1 6 u 3 (23)
where the geometrical boundary conditions, u x (0, t) = 0, ϑ(0, t) = 0 [START_REF] Nayfeh | On the Discretization of Spatially Continuous Systems with Quadratic and Cubic Nonlinearities[END_REF] have been used. Furthermore, the unique non zero strain measure, the curvature χ(s, t) := ϑ , becomes:

χ = u + 1 2 u 2 u ( 25
)
Denoting by t (s, t) and m (s, t) the internal contact force and couple respectively, and b(s, t) and c(s, t) the distributed body forces and couples per unit reference length, respectively, the dynamic balance equations, in Lagrangian form, are expressed as:

t + b = 0, m + x × t + c = 0 (26)
along with the mechanical boundary conditions:

t B = -(k e u B + c e uB ) a y -P d B m B = -c t θB a z (27)
where B denotes evaluation at the end B of the beam and the overdot indicates differentiation with respect to time. Equilibrium in C 0 entails that t 0 (s) = -P a x and m 0 (s) = 0. To obtain the scalar form of the incremental balance equations, the internal contact forces and couples in the current configuration C t are expressed as t(s, t) = [-P + H(s, t)] a x + V (s, t) a y and m(s, t) = M (s, t) a z . Recalling that, for Kirchhoff's beam model, x = d x = cos ϑa x +sin ϑa y , the equations of motions [START_REF] Lacarbonara | Elastodynamics of nonshallow suspended cables: Nonlinear vibration characteristics[END_REF] for the dynamic components, projected onto (a x , a y )-basis, are:

H -mü x = 0, V -mü = 0, M + V cos ϑ + (P -H) sin ϑ = 0 ( 28
)
while the boundary conditions Eq. ( 27) becomes:

H B = P (1 -cos ϑ B ), V B = -(k e u B + c e uB ) -P sin ϑ B , M B = -c t θB (29)

A c c e p t e d m a n u s c r i p t B Eigenvalue Transversality Condition

Let O ≡ (κ 0 , µ 0 ) a point of the parameter space at which the eigenvalue λ is purely imaginary, namely λ 0 := iω 0 (Hopf-point). In order the bifurcation occurs, however, [∂λ/∂ε] = 0, with ε a perturbation of the parameters (not tangent to H). By letting:

λ = λ 0 + ελ 1 + • • • , c = c 0 + εc 1 + • • • (35)
in the eigenvalue problem S(λ; ε)c(ε) = 0, the perturbation equations follow (see also [34]):

ε 0 : S 0 c 0 = 0 ε : S 0 c 1 = S, 0 ε +S, 0 λ λ 1 c 0 ( 36
)
where a comma denotes total differentiation with respect the following variable, and the apex 0 denotes evaluation at C. Enforcing solvability of the ε-order equation, it follows:

λ 1 = - b H 0 S, 0 ε c 0 b H 0 S, 0 λ c 0 (37)
where b 0 is the left eigenvector of S 0 (namely, (S 0 ) H b 0 = 0, with H the transpose conjugate).

Equation (37) represents the sensitivity of λ at O; if [λ 1 ] = 0, then transversality is satisfied.

C Coefficients of the bifurcation equations

The numerical values assumed by the critical parameters, the eigenpairs and the coefficients appearing in the bifurcation equations [START_REF] Luongo | Multiscale analysis of defective multiple-Hopf bifurcations[END_REF] are reported here for the two selected points, HH and HD of Fig. 8.

Double Hopf

The critical values of the parameters are κ 0 = 33.478, µ 0 = 15.602, ξ t = 0.250, ξ e = 0.125.

The frequencies ω i (i = 1, 2) associated with the two interactive modes, and the coefficients p, q, c 1 and c 2 appearing in the relevant (right) critical eigenvectors φ i (s) (Eqs. ( 6) and ( 7))

are given in the Table 1. The values assumed by coefficients c 1 and c 2 in the associated left eigenvectors ψ i (s) (still given by Eq. ( 7)), are also given in the

= (C 1β β + C 1γ γ)A 1 + C 111 A 2 1 Ā1 + C 122 A 1 A 2 Ā2 Ȧ2 = (C 2β β + C 2γ γ)A 2 + C 112 A 1 Ā1 A 2 + C 222 A 2 2 Ā2 (38)
with coefficients reported in the Table 2. In Eq. ( 17),

α iβ = [C iβ ], α iγ = [C iγ ], R ijk = [C ijk ],
while the imaginary parts affect the frequency corrections, not investigated here.

Hopf-divergence

The critical values of the parameters are κ 0 = 36.569, µ 0 = 14.926, ξ t = 0.250, ξ e = 0.125.

The frequencies ω i (i = 1, 2) of the two modes and the parameters p and q are reported in Table 3. Concerning the eigenvectors, φ 1 (s) assumes the following form when indeterminacy has been resolved:

φ 1 (s) = c 1 (cos px -1) + c 2 (sin px -px) (39)
and analogously for ψ 1 (s); Eq. ( 7) still holds for φ 2 (s) and ψ 2 (s). The relevant coefficient c 1 and c 2 are displayed in the Table 3. The complex bifurcation equations are:

Ȧ1 = (C 1β β + C 1γ γ)A 1 + C 111 A 3 1 + C 122 A 1 A 2 Ā2 Ȧ2 = (C 2β β + C 2γ γ)A 2 + C 112 A 2 1 A 2 + C 222 A 2 2 Ā2 (40)
with coefficient given in the Table 4. The same notation adopted for the double-Hopf holds for the coefficients α iβ , α iγ and R ijk of Eq. [START_REF] Luongo | Multiscale analysis of defective multiple-Hopf bifurcations[END_REF].

A c c e p t e d m a n u s c r i p t 12 A c c e p t e d m a n u s c r i p t 1

 121 Linear and Nonlinear Interactions Between Static and Dynamic Bifurcations of Damped Planar Beams Angelo Di Egidio, Angelo Luongo 1 and Achille Paolone Total number of pages: 40 Total number of tables: 4 Total number of figures: Introduction

Figure 2

 2 Figure 2 shows the linear stability diagram of the undamped (ξ e = ξ t = 0) system, where only the (more significant) lower branches of the D-and H-curves have been plotted. The Hopf curve intersects the µ-axis at the well-known µ = 10.025 value of the Beck beam. Then, for increasing κ, the dynamical critical load increases, and the curve tangentially merges at a (degenerate) double zero DZ-point with the upper part of curve D. The disappearance,

Figures 6

 6 Figures6 and 7are now analyzed in more details. In Fig.6a , both damping coefficients was varied, by keeping constant their ratio Ξ := ξ e /ξ t . It was observed, that all the open (small κ) Hopf curves die at the same DZ-point on the divergence curve, as already noted in literature for simple, discrete systems[START_REF] Hermann | Stability of equilibrium of elastic systems subjected to non-conservative forces[END_REF]. The close (larger κ) Hopf curves also intersect the D-curve at different, damping-dependent, HD-points. The influence of the Ξ-parameter is shown in Fig.6b , where κ was kept fixed and Ξ varied. It is observed that the range of ξ e , in which the dynamic bifurcation occurs at a load lower than the static bifurcation, sensibly depends on Ξ. In Figure7a only ξ e was varied, by keeping ξ t fixed. The property of coalescence of DZ-points does not hold any more (as for discrete systems); in addition, for

  a ), for β increasing along path I, the trivial solution is (locally) stable up to point A, where a stable monomodal a 2 -motion takes place. Then, at point B, such motion exchanges stability with a bimodal (a 1 , a 2 )-motion. This disappears at point C, where only an unstable a 1motion, originating at D, exists. Behind D, only the unstable a 2 -motion, originated at A, exists. Figure10 b, relevant to path II, is self-explaining. Figure10c illustrates path III.

A c c e p t e d m a n u s c r i p t

  a ), a divergence of amplitude a 1 occurs at point A. Then, at point B, two bifurcations coalesce: an unstable periodic a 2 -motion and a new, stable, two-components motion (a 1 , a 2 ). Along path II (Fig.12 b), a stable periodic motion around the trivial configuration is triggered C; then, at D, after an unstable static bifurcation, a new equilibrium path arises. Finally, at point E, the motion around the buckled configuration loses stability, and a stable (a 1 , a 2 ) periodic motion around a new equilibrium position takes place.

A c c e p t e d m a n u s c r i p t

  expressing shear-undeformability and inextensibility, respectively (where (•) := d/ds). By taking u := u y as active kinematic variable, from Eqs. (

ACKNOWLEDGMENTS:

  This work was partially supported under the FY 2005-2007 PRIN Grant 'Modelling and experimental tests of the dynamic behavior of flexible structures' from the Italian Ministry of Education, University and Scientific Research.

A c c e p t e d m a n u s c r i p t [ 2 ]

 2 Figure 1: Model.

Figure 2 :

 2 Figure 2: Linear stability diagram for the undamped system.

Figure 3 :Figure 4 :Figure 5 :Figure 6 : 25 Figure 7 :Figure 8 :Figure 9 :

 345625789 Figure 3: Linear stability diagram for the extensionally damped system.

Figure 10 :

 10 Figure 10: Equilibrium paths for the HH bifurcation.

Figure 11 :

 11 Figure 11: Bifurcation diagram and phase-portraits for the HD bifurcation.

Figure 12 :

 12 Figure 12: Equilibrium paths for the HD bifurcation.

A c c e p t e d m a n u s c r i p t

Table 1 .

 1 The complex

	bifurcation equations read:
	Ȧ1
	A c c e p t e d	m a n u s c r i p t

Table 1 :

 1

	φ(s)	ψ(s)

Table 3 :

 3

	φ(s)	ψ(s)

A c c e p t e d m a n u s c r i p t

having neglected the rotatory inertia and expressed the translational inertial forces in the form b = -m(ü x a x + üa y ), where m is the mass density.

In Eqs. ( 28), H and V are reactive internal forces. In order to express the equilibrium equations in terms of the unique configuration variable u, they must be condensed. To this end, from Eq. (28 1 ) the component H is derived, namely:

where the boundary condition Eq. ( 291 ) has been used. To eliminate the V component Eq. ( 283 ) is solved with respect to V , thereby obtaining V = -M / cos ϑ + (H -P ) tan ϑ.

Then, substituting V into Eq. ( 282 ) it follows:

while the boundary conditions Eq. (29 2,3 ) becomes:

A linear constitutive equation is assumed in the form M = EIχ, where EI is the flexural stiffness of the beam. Substituting it in Eq. ( 31), using Eq. ( 22), and Taylor expanding up-to third order, the following equation is obtained: