On the dynamics of Rayleigh beams resting on fractional-order viscoelastic Pasternak foundations subjected to moving loads

Chaos, Solitons and Fractals 93 (2016) 39–47

Contents lists available at ScienceDirect

Chaos, Solitons and Fractals Nonlinear Science, and Nonequilibrium and Complex Phenomena journal homepage: www.elsevier.com/locate/chaos

On the dynamics of Rayleigh beams resting on fractional-order viscoelastic Pasternak foundations subjected to moving loads L.M. Anague Tabejieu, B.R. Nana Nbendjo∗, P. Woafo Laboratory of Modelling and Simulation in Engineering, Biomimetics and Prototypes, Faculty of Science, University of Yaoundé I, P.O. Box 812, Yaoundé, Cameroon

a r t i c l e

i n f o

Article history: Received 30 May 2016 Revised 22 July 2016 Accepted 1 October 2016

Keywords: Fractional-order derivative Rayleigh beams Pasternak foundation Moving loads Horseshoes chaos Beam stability

a b s t r a c t The standard averaging method is used to provide an analytical explanation on the effects of spacing loads, load velocity, order of the fractional viscoelastic property of shear layer material on the amplitude of the beam. The geometric nonlinearity is taken into account in the model. The analysis shows that, when the moving loads are uniformly distributed upon all the length of the structure, it vibrates the least possible. Moreover, as the order of the derivative increases, the resonant amplitude of the beam vibration decreases. In other hand, by means of Melnikov technique, a necessary condition for onset of horseshoes chaos resulting from heteroclinic bifurcation is derived analytically. We point out the critical weight of moving loads and order of the fractional derivative above which the system becomes unstable.

1. Introduction Many structures, such as bridges, runways, rails, roadways, pipelines, etc., can be modelled as a beam structure on a viscoelastic foundation. In the analysis of vibration of beams resting on viscoelastic foundation under a moving load, the beam can be modelled as an Euler–Bernoulli beam [1–4], or as a Rayleigh beam [5] or as a Timoshenko beam [6–8] and the viscoelastic foundation as a Winkler model [9–11] or a Pasternak model [12–14]. The extensive research in this field is summarized in review articles by Fryba [15], Wang et al. [16], and Beskou and Theodorakopoulos [17]. Amongst those investigations, Kargarnovin and Younesian [13,14] used the first order perturbation method to analyse the response of an infinite Timoshenko beam on the viscoelastic foundation under a moving load. To simulate the behaviour of the foundation, Pasternak viscoelastic model was used. This model includes a Kelvin foundation in conjunction with a shear viscous layer [14] or with a shear elastic layer [13]. These studies clearly show that the shear layer material of the foundation can have viscoelastic physical properties. In the context of the bridges behaviour, it is well known that, this layer of material (part ranging between the bridge deck and the foundation) is constituted by some viscoelastic materials such as elastomer. Therefore, the long memory effects of this viscoelas∗

Corresponding author. E-mail addresses: [email protected], [email protected] (B.R. Nana Nbendjo). http://dx.doi.org/10.1016/j.chaos.2016.10.001 0960-0779/© 2016 Elsevier Ltd. All rights reserved.

© 2016 Elsevier Ltd. All rights reserved.

tic materials must be considered through a fractional-order derivative concept. For this aim, we focus in this work on the analytical and numerical analysis of Rayleigh beams subject to uniform moving loads resting on Pasternak foundations considering the shear layer as fractional-order viscoelastic material. Atanackovic et al. [18,19] considered the vibrations of an elastic rod loaded by axial force of constant intensity and positioned on a foundation having fractional order viscoelastic physical properties. The viscoelastic foundation is modelled by using the constitutive equation of Kelvin–Voigt type containing fractional derivatives of real and complex order [18], or by using Zener model of fractional derivatives [19]. The aim of this work is to explore analytically and numerically, the effects of loads number and their spacing, the load velocity, the order of the fractional viscoelastic shear layer material of the Pasternak foundation and its strength on the amplitude of vibration of the beam, and especially the order of the fractional viscoelastic property of the shear layer on the appearance or disappearance of horseshoes chaos through the Melnikov criteria [20– 23]. The paper is organized as follows. Section 2 deals with the derivation of the physical system and their model-building, some mathematical development and the equivalent modal equation of this system. In Section 3, analytical and numerical method (the averaging method [24–27] and Newton–Leipnik methods [28,29] respectively ) are used to analyse the effects of the main parameters of the system on the steady-state amplitude and the stability of the beam system. Section 4 deals with the influence of order

40

L.M. Anague Tabejieu et al. / Chaos, Solitons and Fractals 93 (2016) 39–47

attention is paid to the geometric nonlinearities due to the Euler– Bernoulli law, that states that the bending moment of the beam is proportional to the change in the curvature produced by the load. We call the attention of the readers by taking into account the fact that the geometric nonlinearities due to the axial stiffness and axial (catenary) effects [32] have not been used in this work. To facilitate a compact representation of the governing equation, a window function ε i is employed [39]: εi = 0 when the loads have left the beam and εi = 1 while the loads are crossing the beam. Moreover, QF (x, t) is the foundation-beam interaction force (per unit length of the beam’s axis) which is obtained (including the new term containing fractional derivative) as [13]: Fig. 1. Sketch of a Rayleigh beam on viscoelastic Pasternak foundation under equidistant moving loads. The gravitational forces are represented by arrows P, whose separations are uniform, for the identical speed v0 .

of the fractional viscoelastic shear layer material on the threshold condition for the inhibition of horseshoes chaos in the system. Section 5 is devoted to the conclusion. 2. Beam-foundation model 2.1. Mathematical modelling

θ (x, t ) 

∂ w(x,t ) ∂x

is obtained as:

  2  3 ∂ w(x, t ) ∂ 2 w(x, t ) ∂ 2 ∂ 2 w(x, t ) ρS + EI 2 1− 2 ∂x ∂t2 ∂x ∂ x2 −ρ I



N−1

∂ 4 w(x, t ) + QF (x, t ) = P εi δ [x − xi (t − ti )] 2 2 ∂ x ∂t i=0

(1)

In which S, E, I, ρ , w(x, t) are cross-sectional area of the beam, the modulus of elasticity, cross-sectional moment of inertia, beam material density and the transversal deflection of the beam element, respectively. xi (t − ti ) = v0 (t − id/v0 ) is the position of the ith force at the time t, ti = id/v0 = arriving time of the ith load at the beam, and N is the number of the applied loads. Eq. (1) contains a nonlinear term which may be caused by large curvatures of the beam due to a high deflection [30,31]. So in this work, our

∂ w(x, t ) ∂ 2 w(x, t ) − [μe + μv Dtα ] ∂t ∂ x2

(2)

In which k and c are foundation stiffness and damping coefficients, μe and μv are foundation shear elastic and viscosity coefficients. Dtα is the fractional derivative with order α ∈ (0, 1). There are several definitions for fractional-order derivative, and they are equivalent under some conditions for a wide class of functions. The suitable definition is that of Caputo because it provides initial conditions which physically can be explained [28]. So in Caputo’s sense, the definition of α order derivative of f(t) to t is:

Dtα [ f (t )] =

With reference to Fig. 1, consider a simply supported Rayleigh beam of finite length L, placed on a viscoelastic Pasternak foundation, subjected to a series of lumped loads P with constant interval d moving at the same direction with constant velocity v0 . The loading actions is idealized a real train moving on the bridge. It is assumed that the Pasternak viscoelastic foundation includes a Winkler foundation in conjunction with a shear layer material which is modelled here by using the constitutive equation of Kelvin–Voigt type containing fractional derivative of real order. Fractional hereditary materials involving Kelvin–Voigt units with real order fractional derivatives are analysed in [33,34]. The soil is considered to be infinitely extending beyond the beam and the displacement originating from the beam to propagate along it (see [35–37] ). To take into account the high frequency motion of the beam, a Rayleigh beam correction (up to the second order of the bending angle) [38] is used to refine the theory of Euler– Bernoulli beam for high frequency motion. When the governing equation for the vertical displacement of the beam incorporates the Rayleigh term into the analysis, the correction affects the summing moments produced in the simple Euler–Bernoulli theory. The deformed beam can be described by the transverse deflection w(x, t) and the rotation of the cross section of the beam θ (x, t). Considering the Newton’s second law of motion for an infinitesimal element  of the beam, theequation of motion for the small deformations

QF (x, t ) = kw(x, t ) + c

dα f (t ) 1 = dt α (1 − α )



t 0

f  (u ) du (t − u )α

(3)

where  (z) is Gamma function that satisfies (z + 1 ) = z(z ). Boundary conditions corresponding to the beam shown in Fig. 1, are:

∂ 2 w(x, t ) = 0. ∂ x2 x=0,L

w(x, t )|x=0,L = 0,

(4)

The initial conditions are supposed to be zero when the first moving force enters the beam:

w(x, t )|t=0 = 0,

∂ w(x, t ) =0 ∂ t t=0

(5)

2.2. Modal equations For the analytical purpose, it is convenient to assume an expansion of the transversal deflection w(x, t) in a series form as:

w(x, t ) =

∞

qn (t ) sin

n=1

 nπ x  L

(6)

where qn (t) is the amplitude of the nth mode, and sin (nπ x/L) is the solution of the eigenvalue problem which depends on the boundary conditions of the free oscillations of the beam. After substituting Eqs. (2) and(6) into Eq. (1), multiplying both sides of the resultant equation by the shape function sin (nπ x/L), then integrating with respect to the beam axis x over the length L, and at last considering the following dimensionless variables

χn =

qn , lr

τ = ω0 t

(7)

we obtain for the first mode of vibration [40,41], the dimensionless modal equation (with n = 1 and χ1 = χ (τ )) given by:

χ¨ (τ ) + 2ξ χ˙ (τ ) + 20 χ (τ ) + βχ 3 (τ ) + λDατ χ (τ )

 N−1

dω = P0 εi sin  τ − i 0 v0 i=0 with

(8)

L.M. Anague Tabejieu et al. / Chaos, Solitons and Fractals 93 (2016) 39–47

π v0 2P L3 c L3 , P0 = , ξ = ,  4 Lω0 lr EIπ 2π 2 EIρ [L2 S + Iπ 2 ]  2 3 π lr μv L 2 β=− , λ = ηω0α , η = , 8 L EIπ 2

For that aim, let us assume that the solution of Eq. (12) can be written as

=

20 = 1 + and

ω0 =

π2

(9)

EI

 , ρ L2 S + I π 2

L

χ (τ ) = a(τ ) cos (τ + ϕ (τ )), χ˙ (τ ) = −a(τ ) sin (τ + ϕ (τ )).

k L 2 + μe π 2 kcr L2



L lr = , 2

EIπ 4 kcr = , L4

μecr =

kcr L2

π2

.

(10) According to the Handbook of Mangulis [42], and after considering the configuration of the system where all the applied loads are crossing the beam (εi = 1), the sum of the second member of Eq. (8) becomes N−1



εi sin  τ − i

i=0

d ω0



v0

= sin (τ ) +

 2 sin =

N−1



sin 

i=1

d ω 0 2v0



τ −i

d ω0



v0



sin (N − 1 ) 2dvω0 0



d ω 0 1 − cos v0

 d ω 0 × sin τ − N 2v0







1 a˙ = −  [M1 (a, ϕ ) + M2 (a, ϕ )] sin (τ + ϕ ) 1 aϕ˙ = −  [M1 (a, ϕ ) + M2 (a, ϕ )] cos (τ + ϕ )

M1 (a, ϕ ) = F0N sin (τ ) − G0N cos (τ )



where

− β a3 cos3 (τ + ϕ )

M2 (a, ϕ ) = −λDατ (a cos (τ + ϕ ) )

F0N = P0 G0N =



2P0 sin τ˜0 sin ( (N − 1 )τ˜0 ) sin (N τ˜0 ), 1 − cos (2τ˜0 )

τ˜0 =

dπ . 2L

(13)

From Eq. (12), the dynamics of the beam is that of a Duffing like oscillator with a catastrophic monostable potential since β < 0. This configuration of the potential appears to be more realistic since it may explain best the full dynamics response due to small or high external excitation from a system in engineering science point of view [43]. Noticed that when all the applied loads have left the beam (εi = 0), the second member of Eq. (8) is now equal to zero. Therefore, for the simplicity, we assume in this case that the dynamic response of the system is also equal to zero (χ (τ ) = 0 ) [44]. In the following section, the study of the nonlinear resonance and the stability of the steady-state solution of the considered system will be carried out. 3. Dynamical analysis In this section, a particular attention is focused on the analytical and numerical analysis of the spacing and the velocity of the moving loads, the order of the fractional viscoelastic shear layer material and its strength on the beam response.

(16)

To apply the method [24–27], we average at the period T of which one could select as T = 2π / in the case of periodic function (M1 (a, ϕ )), or T = ∞ in the case of aperiodic one (M2 (a, ϕ )). We obtain the following pair of first order differential equations for the amplitude a(τ ) and the phase ϕ (τ )

aϕ˙ =

λa

 1

 +

2 sin τ˜0 sin ( (N − 1 )τ˜0 ) 1+ cos (N τ˜0 ) , 1 − cos (2τ˜0 )



+ a 2 − 20 cos (τ + ϕ ) + 2ξ a sin (τ + ϕ )

So Eq. (8) can be rewritten in the alternative form as

(12)

(15)

where

a˙ = −ξ a −

(11)

(14)

where the amplitude a(τ ) and the phase ϕ (τ ) are slow-varying functions of time τ . Substituting Eq. (14) into Eq. (12), we obtain

and

χ¨ (τ ) + 2ξ χ˙ (τ ) + 20 χ (τ ) + βχ 3 (τ ) + λDατ χ (τ ) = F0N sin (τ ) − G0N cos (τ )

41

−

2

α−1 sin



a 2 − 20 2

 απ  2



+

3β 3 + a 8

1 [G0N sin ϕ − F0N cos ϕ ] 2

+

λa 2

α−1 cos

(17)

 απ  2

1 [G0N cos ϕ + F0N sin ϕ ] 2

(18)

Now, we study the steady-state solution, which is more important and meaningful in vibration engineering. By putting a = A0 , ϕ = 0 and a˙ = 0, ϕ˙ = 0, and After eliminating 0 from Eq. (18), we find that the amplitude A0 of the oscillatory state satisfies the following nonlinear equation:

  9 2 6 3 β A0 − β 1 (α )A40 + 21 (α ) + 22 (α ) A20 = F02N + G20N 16 2

(19)

with

 απ    1 (α ) = 2 − 20 − λα cos 2  απ  α 2 (α ) = 2ξ + λ sin

(20)

2

This equation has more than one steady-state solution for some parameters. An interesting observation is the dependence of the oscillations amplitude upon the beam parameters (natural frequency 20 and nonlinear component β ), the loads traffic (loads weights intensity P0 , loads number N, spacing loads d and driving frequency  ) and of the foundation parameters (damping ξ , the order of the fractional derivative α and its strength λ). Next, we study the stability of the steady-state solution by using the method of Andronov and Witt [45]. Let a = A0 + a , ϕ = 0 + ϕ and substituting them into Eqs. (17) and (18) yields



3β 2 d a 2 ( α ) A a + 0  1 ( α ) − =− A ϕ dτ 2 2 4 0 and

(21)



3.1. Approximate analytical solution and its stability analysis

d ϕ 1 9β 2 2 ( α ) = A −  1 ( α ) a − ϕ dτ 2 A 0 4 0 2

We start by using the averaging method [24–27], which provides an analytical approximate solution and thus permits to detect the effects of the main parameters on the system response.

where 1 (α ) and 2 (α ) are given by Eq. (20). The stability of the steady-state solution is determined by the eigenvalues of the corresponding Jacobian matrix of Eqs. (21) and (22). The corresponding

(22)

42

L.M. Anague Tabejieu et al. / Chaos, Solitons and Fractals 93 (2016) 39–47 Table 1 Properties of the beam, foundation and moving load. Item Beam Length Young’s modulus (steel) Cross-sectional area Mass density Moment of inertia

Notation

Value

L E S

628.1 m 200 MPa 4.8 m2 7850 kg /m3 12 m4

ρ

I

Foundation Stiffness Viscous damping shear viscosity coefficient Shear stiffness coefficient Moving load Load Mean velocity

μe

202.66 N/m2 192 Ns/m2 5 × 106 Ns 5 × 105 N

P v0

350 KN 30 m/s

k c

μv

Fig. 2. Amplitude response of the system A0 as function of the driven frequency  for different values of the fractional-order α . N = 16, d = L/(N − 1 ).

eigenvalues  are the roots of

 2

 2 ( α ) 2 ( α ) 1 3β 2  + + A − 1 ( α ) +  2 4 4 0

 9β 2 × A 0 − 1 ( α ) = 0 2

(23)

4

Since 0 < α ≤ 1, then 2 (α ) > 0 and the instability condition for the steady-state solution is found by using the Routh-Hurwitz criterion [46] as



2 ( α ) 2

2

+



1 3β 2 A − 1 ( α ) × 4 4 0



9β 2 A − 1 ( α ) < 0 4 0 (24)

This condition keeps the real parts of the eigenvalues positive. The asymptotically stable solution may occur in the opposite case. 3.2. Numerical analysis The physical and geometrical properties of the beam are listed in Table 1 [47,48]. Thus, the dimensionless parameters of Eq. (12) are:

ξ = 0.012; β = −0.926; P0 = 0.002;  = 0.751; λ = 0.832 × (0.2)α ; 20 = 1.143 In order to verify the precision of the analytical solutions, we first solve numerically Eq. (12) using the Newton–Leipnik algorithm [28,29] by considering the following definition of the fractional order derivative:

Dατ [χ (τn )] ≈ h−α

n

C αj χ



τn − j



(25)

j=0

where h is the integration step and the coefficients C αj satisfy the following recursive relations:

C0α = 1,



C αj = 1 −



1+α α C j−1 j

(26)

Second, we display in some figures the effects of the main parameters of the proposed model. For example, Fig. 2 shows the effects of the order of the derivative on the amplitude of vibration of the beam. This graph also shows a comparison between the results from the mathematical analysis (curve with dotted line) and the results obtained from numerical simulation of Eq. (12) (curve with circles for α = 0.25) using the definition of Eq. (25). The match between the results shows a good level of precision of the approximation made in obtaining Eq. (19). This figure also reveals that as the order of the derivative increases, the resonant amplitude of the

beam vibration decreases. Similar results were obtained by Shen et al. [26] who analyzed the periodic vibration of a Duffing oscillator with additional fractional damping and showed that increasing the order of the derivative (damping term) results in a decrease of the amplitude of vibration. For small order derivative, the bending degree of the curves becomes severe, which leads to multiple amplitude solutions and the appearance of the unstable solution. In Fig. 3, we have plotted the evolution of the amplitude of vibration of the beam A0 as a function of the shear viscosity coefficient η for several orders of the derivative. It clearly shows that the system is more stable for the highest order of the derivative. The multivalued solution appears for small order and disappears progressively as the order increases. Fig. 4 shows first how an increase in the number of the moving loads affects the amplitude of vibration of the beam. It is observed that as the value of N increases, the amplitude of vibration at the resonant state merely increases. Second, the effect of spacing loads is investigated. It is found that when the moving loads are uniformly distributed upon all the length of the structure, it vibrates the least possible. The so-called force-response curve is depicted in Fig. 5. Here, multiple and up to three coexisting solutions (the dotted and circles lines correspond to stable and unstable branches, respectively) are observed. In particular, there are three coexisting solutions for P01 < P0 < P02 , and exactly one solution branch outside this region of bistability. The stable (dotted line) and unstable (circles line) branches merge at P0 = P01 and P0 = P02 . For the practical convenience, it is important to note that, the bistability region (hatched zone of Fig. 5) gives the range of dangerous weight of the moving loads for the beam safety. In this section, we have studied the effects of some main parameters such as α , η, N and d on the periodic solution of the system. The following section, will investigate the peculiar effects of fractional order of the derivative and the shear viscosity coefficient on the possible appearance of the horseshoes chaos by using the Melnikov theory. 4. Horseshoes chaos caused or controlled by the order of the derivative: Melnikov analysis In the present section, we apply the Melnikov method [20– 23] to detect analytically the effects of the fractional order of the derivative and the shear viscosity coefficient on the threshold condition for the inhibition of smale horseshoes chaos in the system and on the fractal basin boundaries. Many researchers have begun to investigate the chaotic dynamics of fractional-order systems

L.M. Anague Tabejieu et al. / Chaos, Solitons and Fractals 93 (2016) 39–47

Fig. 3. The steady-state amplitude of the beam A0 as function of the shear viscosity coefficient η for several order of the derivative. N = 16,  = 1.14, d = L/(N − 1 ).

Fig. 4. Vibration amplitude of the beam A0 for different values of the loads number N versus the dimensionless spacing loads d/L. α = 0.8,  = 1.14.

43

44

L.M. Anague Tabejieu et al. / Chaos, Solitons and Fractals 93 (2016) 39–47

Fig. 5. Force-response curve for α = 0.8,  = 0.9, N = 3, d = L/(N − 1 ). Fig. 6. Phase space trajectories.

[49,50]. A peculiar attention is put on the work of Oumbé et al. [49] who, based to the Melnikov method, demonstrated that the order and strength of the fractional viscoelastic property of the flexible material can be effectively used to control chaos in a system. Our main objective in this section is to use the same method to demonstrate that the order of the fractional viscoelastic shear layer material can be used either to control chaos in the beam system, or to cause chaos. To perform the Melnikov analysis, we introduce the small parameter ε into the differential Eq. (12) and rewrite the governing system as

dU = F [U ] + ε G[U, τ ] dτ



G=



χ

y = χ˙

,

F=

−2ξ y − λ ∗

∗

y −20 χ − βχ 3





0

(28)

Dατ χ + F0∗N sin τ − G∗0N cos τ

1 1 β H (χ , y ) = y2 + 20 χ 2 + 2 2 4

χ4

(29)

Since β < 0, the system has three equilibrium points that area center point χes = (0, 0 ) and two saddle points χeu± =

20 β , 0 ). Since phase trajectories cannot cross the center, the

unperturbed system shows a stable periodic motion. These trajectories are shown in Fig. 6. The saddle points surrounding the center are connected by heteroclinic orbits that satisfy the following equation:



χhet (τ ) = ±0 . − β1 tanh χ˙ het (τ ) = yhet (τ ) = ±

2 0.





0 √ . 2

τ

− 21β



sec h

2



0 √

2

.τ





+∞ −∞

= −2ξ

F [Uhet (τ )] ∧ G[Uhet (τ ), τ + τ0 ]dτ



+ F0N

+∞

y2het (τ ) dτ − λ

−∞  +∞

− G0N

−∞

+∞

−∞

yhet (τ )Dατ [χhet (τ )] dτ

yhet (τ ) sin (τ + τ0 )dτ

−∞ +∞





yhet (τ ) cos (τ + τ0 )dτ

(31)

When the Melnikov function has simple zero point, the stable manifold and the unstable manifold intersect transversally, chaos in the sense of Smale horseshoe transform occurs. So let MD (τ0 ) = 0, one concludes that horsehoes chaos appears when

with ξ = ε ξ ∗ , λ = ε λ∗ , F0N = ε F0∗N , G0N = ε G∗0N and ε the perturbation parameter. In the rest of the paper, the stars (∗) are removed for simplicity. For ε = 0, the system of Eq. (27) is the Hamiltonian system with Hamiltonian function

(± −

M D ( τ0 ) =

(27)

where the vector fields U, F and G are given by

U=

follows

(30)

The Melnikov theory defines the condition for the appearance of the so-called transverse intersection points between the perturbed and the unperturbed separatrix or the appearance of the fractality on the basin of attraction. This theory can be applied in the case of Eq. (27) by using formula given by Wiggins [20] as

P0 ≥ P0 cr

√ 2ξ 0 I1 − λ 2Iα = 

 −2β 2 sin τ ˜ sin ( N − 1 ) τ ˜ ( ) 0 0 .(cos(N τ˜0 ) + sin(N τ˜0 )) 0 I2 1 + 1 − cos (2τ˜0 ) (32)

where

I1 =

√ 4 2 3 0

2π 

 π √ 0 2  

 +∞ 0 0 Iα = sec h2 √ τ Dατ tanh √ τ dτ I2 =

20 . sinh

−∞

2

2

(33)

The integral Iα is evaluated by using the numerical formulas for the right and left side of the fractional derivative [49]. Here P0 cr is the threshold amplitude for the onset of horseshoes chaos in the system. This threshold condition is plotted in Fig. 7 as a function of the driving frequency  for different values of the fractional order α (Fig. 7(a) and(b)), as function of the viscosity coefficient η (Fig. 7(c)) and as function of the fractional order (Fig. 7(d)). We observe that the area below the curves (dotted line) indicates the region for regular solutions while above them the corresponding solutions have sensitivity to initial conditions, chaotic transient and fractal basin boundaries. The graph also shows that increasing the order of derivative contributes first to lower the threshold

L.M. Anague Tabejieu et al. / Chaos, Solitons and Fractals 93 (2016) 39–47

45

Fig. 7. Critical weight P0 cr for the appearance or disappearance of horseshoes chaos as a function of the system parameters. N = 6, d = L/(N − 1 ).

boundary for onset of chaos and then enlarges the possible chaotic domain (Fig. 7(a)), and second to rise it (reduction of the chaotic domain) (Fig. 7(b)). This observation is well illustrated in Fig. 7(d) where we demonstrate that for 0 < α < 0.48, the threshold decreases and for 0.48 < α < 1 , it increases. The same investigations are made on the case of Fig. 7(c). To validate the accuracy of the proposed analytical predictions, we solve numerically Eq. (12) using Newton–Leipnik method [28,29] to display the shape of the basin of attraction. A particular characteristic of the Melnikov chaos is the fractality of the basin of attraction and the resulting unpredictability due to the dependence on the initial conditions. Moreover, Awrejcewicz et al. [51] proved a dependence on the fractal structure of a basin of attraction on the occurrence of heteroclinic or homoclinic bifurcation taking as an example the Duffing oscillator. To check the effects of the fractional order α on the performance of our system, Fig. 8 is plotted and the metamorphoses of the basin of attraction around the central equilibrium point (0, 0) are observed. It is clear that the fractal structure of the basin appears for the lowest (α = 0.1 ) and highest (α = 0.86 ) values of the fractional order. While the regular shape is observed for the intermediate values. These observations had already been predicted by the previously analytical developments. In this section, we have studied the effects of fractional order of the derivative and the shear viscosity coefficient on the possible appearance or disappearance of the horseshoes chaos in the beam system. The results obtained here give more precision about the effect of the fractional order on the beam stability investigated previously in Section 3. It is shown how the extreme values (lowest and highest values) of the fractional order have negative effects on the system stability and positive for the intermediate values. Thus, proper selection of this fractional order can be contributed to con-

trol chaos in a system or to cause chaos. These results also explain well why it was necessary to introduce the fractional nature of the shear layer material in the model of a viscoelastic Pasternak foundation. 5. Conclusion We have presented an analytical and numerical solution for the dynamic response of a Rayleigh beam resting on viscoelastic Pasternak foundations considering the shear layer as viscoelastic material having fractional-order property and subjected to the passing of series of equidistant moving loads with constant velocity. The proposed system has been modelled by the full partial differential equation that we have thereafter reduced to a single fractional order nonlinear one-dimensional equation. Based on averaging method, the approximate analytical solution has been obtained and we have studied the effects of some main parameters of the system such as the fractional-order, the shear viscosity parameter of the shear layer foundation, the loads number and their spacing on the vibration amplitude of the beam and on its stability. The analysis leads us to the conclusion that, as the order of the derivative increases, the resonant amplitude of the beam vibration decreases and, when the moving loads are uniform distributed upon all the length of the structure, it vibrates the least possible. Finally, The comparisons of the approximate analytical solution with the numerical one verify the correctness and satisfactory precision of the presented results. The analysis has also allowed to estimate the condition for the possible appearance of Melnikov chaos in the system. The main conclusion is that this condition depends strongly on the order of the derivative. For example, we have demonstrated that increas-

46

L.M. Anague Tabejieu et al. / Chaos, Solitons and Fractals 93 (2016) 39–47

Fig. 8. Basins of attraction showing the confirmation of the analytical prediction for N = 6,  = 1.14, P0 = 0.002, d = L/(N − 1 ).

ing the order of derivative contributes first to lower the threshold boundary for the onset of chaos and then enlarges the possible chaotic domain, and second to rise the threshold (reduction of the chaotic domain). More precisely, proper selection of the shear layer of the foundation that having fractional order viscoelastic material can be contributed to suppression of chaos in a system. The analytical results have been validated numerically via the basin of attraction. Acknowledgements Part of this work was completed during a research visit of Prof Nana Nbendjo at the University of Kassel in Germany. He is grateful to the Alexander von Humboldt Foundation for financial support within the Georg Forster Fellowship. References [1] Ding H, Chen QL, Yang SP. Convergence of galerkin truncation for dynamic response of finite beams on nonlinear foundations under a moving load. J Sound Vib 2012;331:2426–42. [2] Sun L, Luo F. Steady-state dynamic response of a bernoulli-euler beam on a viscoelastic foundation subject to a platoon of moving dynamic loads. J Vib Acoust 2008;130. 051002–1 [3] Dimitrovová Z. A general procedure for the dynamic analysis of finite and infinite beams on piece-wise homogeneous foundation under moving loads. J Sound Vib 2010;329:2635–53. [4] Pinto OC, Gonçales PB. Active non-linear control of buckling and vibrations of a flexible buckled beam. Chaos Solitons Fract 2010;14:227–39. [5] Hryniewicz Z. Dynamics of rayleigh beam on nonlinear foundation due to moving load using adomian decomposition and coiflet expansion. Soil Dyn Earthq Eng 2011;31:1123–31. [6] Kargarnovin MH, Younesian D, Thompson DJ, Jones CJC. Response of beams on nonlinear viscoelastic foundations to harmonic moving loads. Comput Struct 2005;83:1865–77.

[7] Sapountzakis EJ, Kampitsis AE. Nonlinear response of shear deformable beams on tensionless nonlinear viscoelastic foundation under moving loads. J Sound Vib 2011;330:5410–26. [8] Chen YH. Response of an infinite timoshenko beam on a viscoelastic foundation to a harmonic moving load. J Sound Vib 2001;241:809–24. [9] Mathews PM. Vibrations of a beam on elastic foundation. J Appl Math Mech 1958;38:105–15. [10] Haitao Y, Yuan Y. Analytical solution for an infinite euler-bernoulli beam on a visco-elastic foundation subjected to arbitrary dynamic loads. J Eng Mech 2014;140(3):542–51. [11] Achenbach JD, Sun C. Moving load on a flexibly supported timoshenko beam. Int J Solids Struct 1965;1:353–70. [12] Yang Y, Ding H, Chen L-Q. Dynamic response to a moving load of a timoshenko beam resting on a nonlinear viscoelastic foundation. Acta Mech Sin 2013;29(5):718–27. [13] Younesian D, Kargarnovin MH. Response of the beams on random pasternak foundations subjected to harmonic moving loads. J Mech Sci Tech 2009;23:3013–23. [14] Kargarnovin MH, Younesian D. Dynamics of timoshenko beams on pasternak foundation under moving load. Mech Res Commun 2004;31:713–23. [15] Fryba L. Vibration of solids and structures under moving loads. London: Thomas Telford; 1999. [16] Wang YH, Tham LG, Cheung YK. Beams and plates on elastic foundations: a review. Progr Struct Eng Mater 2005;7(4):174–82. [17] Beskou ND, Theodorakopoulos DD. Dynamic effects of moving loads on road pavements: a review. Soil Dyn Earthq Eng 2011;31:547–67. [18] Atanackovic TM, Janev M, Pilipovic S, Zorica D. Vibrations of an elastic rod on a viscoelastic foundation of complex fractional kelvin-voigt type. Meccanica 2015;50(7):1679–92. [19] Atanackovic TM, Stankovic B. Stability of an elastic rod on a fractional derivative type of foundation. J Sound Vib 2004;277:149–61. [20] Wiggins S. Introduction to applied Nonlinear dynamical systems and chaos. New York: Springer; 1990. [21] Melnikov VK. On the stability of the center for some periodic perturbations. Trans Moscow Math Soc 1963;12:1–57. [22] Nana Nbendjo BR, Woafo P. Active control with delay of horseshoes chaos using piezoelectric absorber on a buckled beam under parametric excitation. Chaos Solitons Fract 2007;32:73–9. [23] Fengxian A, Fanggi C. Bifurcations and chaos of the nonlinear viscoelastic plates subjected to subsonic flow and external loads. Chaos Solitons Fract 2016;91:78–85.

L.M. Anague Tabejieu et al. / Chaos, Solitons and Fractals 93 (2016) 39–47 [24] Hayashi C. Nonlinear oscillations in physical sSystems. New York: McGraw-Hill; 1964. [25] Sanders JA, Verhulst F, Murdock J. Averaging methods in nonlinear dynamical systems. New York: Springer Science & Business Media; 2007. [26] Shen Y, Yang Xing H, Ma H. Primary resonance of duffing oscillator with two kinds of fractional-order derivatives. Int J Non-Linear Mech 2012;47:975–83. [27] Deng J, Xie W-C, Pandey MD. Stochastic stability of a fractional viscoelastic column under bounded noise excitation. J Sound Vib 2014;333:1629–43. [28] Debnath L. Recent applications of fractional calculus to science and engineering. Int J Math Math Sci 2003;54:3413–42. [29] Petráš I. Fractional-order nonlinear systems: modeling, analysis and simulation. Beijing: Higher Education Press; 2011. [30] Oumbé Tékam GT, Tchawou Tchuisseu EB, Kitio Kwuimy CA, Woafo P. Analysis of an electromechanical energy harvester system with geometric and ferroresonant nonlinearities. Nonlinear Dyn 2014;76(2):1561–8. [31] Kovacic I, Brennan MJ. The duffing equation: nonlinear osscillators and their behaviour. United Kingdom: John Wiley & Sons; 2011. [32] Nayfeh AH, Mook DT. Non linear oscillations. New York: Wiley-Interscience; 1979. [33] Paola MD, Pinnola FP, Zingales M. A discrete mechanical model of fractional hereditary materials. Meccanica 2013;48(7):1573–86. [34] Eldred LB, Baker WP, Palazotto AN. Kelvin-voigt vs fractional derivative model as constitutive relations for viscoelastic materials. AIAA J 1995;33:547–50. [35] Nobili A. Variational approach to beams resting on two-parameter tensionless elastic foundations. J Appl Mech 2012;79(2):021010. [36] Nobili A, Radi E, Lanzoni L. A cracked infinite kirchhoff plate supported by a two parameter elastic foundation. J Eur Ceram Soc 2014;34(11):2737–44. [37] Lanzoni L, Nobili A, Radi E, Sorzia A. Axisymmetric loading of an elastic– plastic plate on a general two-parameter foundation. J Mech Mater Struct 2015;10(4):459–79. [38] Han SM, Benaroya H, Wei T. Dynamics of transversely vibrating beams using four engineering theories. J Sound Vib 1999;225(5):935–88.

47

[39] Nikkhoo A, Hassanabadi ME, Azam SE, Amiri JV. Vibration of a thin rectangular plate subjected to series of moving inertial loads. Mech Res Commun 2014;55:105–13. [40] Kitio Kwuimy CA, Nana Nbendjo BR, Woafo P. Optimization of electromechanical control of beam dynamics: Analytical method and finite differences simulation. J Sound Vib 2006;298:180–93. [41] Liu D, Xu W, Xu Y. Stochastic response of an axially moving viscoelastic beam with fractional order constitutive relation and random excitations. Acta Mech Sin 2013;29(3):443–51. [42] Mangulis V. Handbook of series for scientists and engineers. New York: Academic Press; 1965. [43] Nana Nbendjo BR, Woafo P. Modelling of the dynamics of euler’s beam by φ 6 potential. Mech Res Commun 2011;38:542–5. ´ [44] Sniady P, Biernat S, Sieniawska R, Zukowski S. Vibrations of the beam due to a load moving with stochastic velocity. Prob Eng Mech 2001;16(1):53–9. [45] Andronov AA, Witt A. Towards mathematical theory of capture. Archiv fur Electrotechnik 1930;24:99–110. [46] Awrejcewicz J, Krysko VA. Introduction to asymptotic methods. New York: Chapmanand Hall, CRC Press, Boca Raton; 2006. [47] Karoumi R. Some modeling aspects in the nonlinear finite element analysis of cable supported bridges. Comput Struct 1999;71:397–412. [48] Yau JD, Yang YB. Vibration reduction for cable-stayed bridges travelled by high-speed trains. Finite Elem Anal Des 2004;40:341–59. [49] Oumbé Tékam GT, Kitio Kwuimy CA, Woafo P. Analysis of tristable energy harvesting system having fractional order viscoelastic material. Chaos 2015;25:013112. [50] Takougang Kingni S, Nana B, Ngueuteu Mbouna GS, Woafo P, Danckaert J. Bursting oscillations in a 3d system with asymmetrically distributed equilibria: mechanism, electronic implementation and fractional derivation effect. Chaos Solitons Fract 2015;71:29–40. [51] Awrejcewicz J, Holiske M. Smooth and non-smooth high dimensional chaos and the Melnikov type methods. Singapore: World Scientific; 2007.