Free Vibrations with Large Amplitude of Axially Loaded Beams on an Elastic Foundation Using the Adomian Modified Decomposition Method

Analytical solutions describing free transverse vibrations with large amplitude of axially loaded Euler–Bernoulli beams for various end restrains resting on a Winkler one-parameter foundation are obtained using the Adomian modified decomposition method (AMDM).)e AMDM allows the governing equation to become a recursive algebraic equation, and, after some additional simple mathematical operations, the equations can be cast as an eigenvector problem whose solution results in the calculation of natural frequencies and corresponding closed-form series solution of the mode shapes. Important to the use of the Adomian modified decomposition method is the treatment of the nonlinear Fredholm integral coefficient, which forms part of the governing equation. In addition to the calculation of natural frequencies and mode shapes, investigations are made of the effects on the free vibrations of the Winkler parameter and of increasing the axial loading.


Introduction
Uniform slender beams resting on an elastic foundation, while subjected to axial loading, are common in structural systems undergoing actual operating conditions.Analysis of such systems, both linear and nonlinear, has been of interest to civil and railway engineering.For example, when ambient temperatures increase, rails and concrete slabs, often used in urban transport systems, tend to expand, so causing compressive in-plane forces, leading to changes in natural frequencies and eventually to buckling.In-plane compressive forces are also found in prestressed beams.If the amplitudes of the vibrations remain small, the governing equation is usually in the form of a linear differential equation which is relatively simple to solve.However for large amplitude vibration, nonlinear terms are introduced into the governing equation which needs to be treated.Rails and concrete slabs often rest on foundations generally classified as elastic, viscoelastic, Winkler, and Pasternak.
Boundary value problems (BVPs) have been the subject of several analytical methods, recently developed, to calculate beams with relatively simple configurations.e variational iteration method (VIM) often attributed to He [1] is a modification of a general Lagrange multiplier method and has been used as a powerful tool for solving ordinary differential equations [2,3].Another recent method developed is the homotopy perturbation method (HPM) [4,5] which has been used for problems involving nonlinear differential equations.Less recent methods used to investigate the vibration problem for nonuniform Euler-Bernoulli beams have been the Rayleigh-Ritz method [6], closed-form solutions [7], and Green's function method [8].Several methods have used the Frobenius series [9] and also by discretizing the beam into beam elements [10].ere have been some early studies of vibrating beams under axial loading [11][12][13], where the effect of increasing the axial loading on the mode shapes and natural frequencies of the beam was investigated.ere has been some work already done for large amplitude vibration.Bhashyam and Prathap [14] used the Galerkin finite element method to study nonlinear vibration, and Özkaya [15] calculated the response of a beam mass system with clamped ends by applying a method known as the method of multiple scales.
For beams resting on foundations, an understanding of the beam-foundation interaction is needed.e foundation increases resistance to movement and can significantly change the modal characteristics of the beam.Many practical cases in engineering related to foundation-beam interaction can be modelled by assuming the beam resting on an elastic foundation with the Winkler elastic foundation model [16] used extensively.
is model assumes the foundation to be made up of an infinitely many closedspaced linear springs and is a one-parameter model.A limitation of this model is the assumption that there is no interaction between the springs.To overcome this defect, several two-parameter models have been suggested, such as Filonenko-Borodich, Pasternak, and Vlasov and Leontiev foundation models [17].Studies using a constant Winkler foundation can be found in the literature [18][19][20].
In the present work, the adomian modified decomposition method [21,22] is utilized to calculate free transverse vibration characteristics of axially loaded Euler-Bernoulli beams with various end restrains, resting on a Winkler one-parameter foundation.e method is chosen as it has proved efficient and accurate [23,24] for solving linear and nonlinear differential equations, and it has the advantage of computational simplicity.In addition, it does not involve linearization, discretization, perturbation, or a priori assumptions, which may alter the physics of the problem considered [21].For the AMDM, the solution is considered to be the sum of an infinite series with rapid convergence [25].Using the AMDM, the governing differential equation becomes a recursive algebraic equation and the boundary conditions become simple algebraic frequency equations, which are suitable for symbolic computation.After some simple algebraic operations on the frequency equations for any i th natural frequency, the closed-form series solution of any i th mode shape can be obtained.Calculations are made for clamped-free and clamped-clamped boundary conditions together with an investigation of the effects of increasing the axial loading and Winkler parameter on the natural vibrations.

Principle of Adomian Modified Decomposition Method (AMDM)
e basic theory of AMDM is briefly stated here.Consider the equation where F represents a general nonlinear ordinary differential operator involving both linear and nonlinear parts and g(x) is a given function.e linear terms in Fy are decomposed into Ly + Ry, where L is an invertible operator, which for AMDM is taken as the highest-order derivative, and R is the remainder of the linear operator.Equation ( 1) can now be written as where Ny represents the nonlinear terms of Fy and equation ( 2) represents an initial value or boundary value problem.
On solving for Ly, equation ( 2) can be transformed to where Φ is the integration constant and LΦ � 0 is satisfied.
To use the AMDM, y is decomposed into the infinite sum of a convergent series: and the nonlinear term is decomposed using Adomian polynomials, A m : e function g(x) can also be decomposed as Putting equations ( 4)-( 6) into equation (3) gives e coefficient C m can be calculated using a recurrence relation and the power series solutions of linear homogeneous differential equations in initial value problems yield simple recurrence relations for the coefficient C m .In practice, the coefficients cannot be determined exactly, and the solutions can only be approximated by a truncated series

Mathematical Formulation
In this work, a uniform beam under axial load while resting on a Winkler foundation is considered.As shown in Figure 1, the beam has a length, l, a uniform rectangular cross section, A, a cross-sectional moment of inertia, I, and the beam is considered as made of isotropic material with a modulus of elasticity, E, and density, ρ. e model for the foundation is the relatively simple Winkler model whose stiffness changes along the beam length and is a function of the spatial coordinate along the beam in the x direction.
According to the theory of structural vibrations [26,27], on using the Euler-Bernoulli beam model, the strain energy induced by a large displacement amplitude is given by 2 Shock and Vibration e large amplitude of the vibrations necessitates the inclusion of the nonlinear term shown in equation (8).
Here, u and w are the axial and transverse displacements, respectively, and k w (x) is the foundation stiffness coefficient.e kinetic energy is given by e external work done by the axial load is By invoking Hamilton's principle and using the Lagrangian of the system, On substituting equations ( 8)-( 10) into equation ( 11), the following governing equation can be obtained after eliminating axial displacement: where k 0 is a constant and g(x) is a function of the spatial coordinate along the beam length.
According to modal analysis for harmonic-free vibration, w(x, t) can be separable in space and time as where ϕ(x) is the modal deflection and h(t) is a harmonic function of time t.If ω denoted the circular frequency of and the eigenvalue problem of equation ( 12) reduces to Equation ( 14) is now made nondimensional using and becomes

Boundary Conditions
Two cases are considered in this work, namely, beams which are clamped-clamped and clamped-free, respectively: For the clamped-clamped case, the boundary conditions at X � 0 and X � 1 are For the clamped-free case, the boundary conditions at X � 0 and X � 1 are It is convenient for the AMDM to describe boundary conditions in terms of rotational and translational flexible ends as shown in Figure 2.
e boundary conditions can be written in the dimensionless form as where the coefficients are nondimensionalized as (20)

Application of the Adomian Modified Decomposition Method (AMDM)
According to the AMDM, ϕ(X) in equation ( 16) can be expressed as an infinite series, i.e., where the unknown coefficients, C m , are determined recurrently.If a linear operator G ≡ d 4 /dX 4 is used, then the inverse operator of G is a four fold operator defined as Equation ( 16) now can be written as where Φ(X) is a polynomial depending on the boundary conditions.

Nonlinear Fredholm Integral
Coefficient.e terms in equation ( 23) are linear except for the nonlinear Fredholm integral coefficient shown in equation (24).Before continuing with the main solution method, the nonlinear term will be treated first through the use of appropriate Cauchy products.Consider Use is now made of the term Differentiating On setting X 0 � 0, So As X 0 � 0 in the current analysis, then Combining the two strands of analysis gives where 5.2.Linear and Nonlinear Terms Combined.Using the analysis of the last subsection and ϕ(X) where 6 is the initial term.e recurrence relation for the coefficients C m can now be stated as and for m ≥ 4 as e coefficients C m can be found from the recurrence equations (34) and (35), and the solution for ϕ(X) is calculated using equation (33).
e series solution is ϕ(X) �  ∞ m�0 C m X m , although all of the coefficients C m cannot be determined, and thus, the solutions must be approximated by the truncated series  n−1 m�0 C m X m and successive approximations are ϕ [n] (X) �  n−1 m�0 C m X m , as n increases and the boundary conditions are met.us, ϕ [1] (X) � C 0 , ϕ [2] (X) � ϕ [1] (X) + C 1 X, ϕ [3] (X) � ϕ [2] (X) + C 2 X 2 , and ϕ [4] (X) � ϕ 3 (X) + C 3 X 3 serve as approximate solutions with increasing accuracy as n ⟶ ∞. e four coefficients C j (j � 0, 1, 2, 3) depend on the boundary conditions used (either equation ( 17) or (18)).For example, for the clamped-free boundary conditions at X � 0, the two coefficients C 0 and C 1 can be chosen as arbitrary constants, and the other two coefficients C 2 and C 3 can be expressed as functions of C 0 and C 1 .us, from equations ( 19) and (34), the following is obtained: us the initial term Φ(X) is only a function of C 0 and C 1 , and from the recurrence relation of equation ( 36), the coefficients C m (m ≥ 4) are functions of C 0 , C 1 , and λ.By substituting ϕ [n] (X) into the boundary conditions of equation ( 19) when X � 1, we have For nontrivial solutions of C 0 and C 1 , the frequency equation is given as e i th estimated eigenvalue λ [n]  (i) corresponding to m is obtained from equation (38), i.e., the i th estimated dimensionless natural frequency Ω [n]  n(i) � ��� λ [n]   (i) is also obtained and n is determined by where is the i th estimated dimensionless natural frequency corresponding to n − 1 and ε is a preset sufficiently small value.If equation (39) is satisfied, then Ω [n]  n(i) is the i th dimensionless natural frequency Ω n(i) .By substituting Ω [n]  n(i) into equation (37), and all of the other coefficients C m can be obtained from equations ( 34) and (35).Furthermore, the i th mode shape ϕ [n]   i corresponding to the i th eigenvalue Ω [n]  n(i) is obtained by where in which λ is substituted by λ i and ϕ [n]   i is the i th eigenfunction corresponding to the i th eigenvalue λ i .By normalizing equation (41), the i th normalised eigenfunction is defined as Shock and Vibration where ϕ [n] i (X) is the i th mode shape function of the beam corresponding to the i th natural frequency ω [n]  i � . is general theory is now applied to a uniform Euler-Bernoulli beam under different boundary conditions.

Clamped-Free Uniform Beam.
e first case considered is the clamped-free uniform beam resting on an elastic foundation and experiencing axial compressive force e case was chosen to test (validate) the accuracy of the AMDM.
e boundary conditions are as given in equations ( 18) and ( 19) with the spring constants becoming When X � 0, the first two boundary conditions of equation ( 19) yield the relationships shown in equation ( 36), and when X � 1 by substituting ϕ [n] (X) �  n−1 m�0 C m X m into the last two boundary conditions of equation ( 19), the following two algebraic equations (written in full) involving C 0 and C 1 are obtained: e case of the clamped-free uniform beam without an elastic foundation support or without any axial force and with small vibration amplitude was first chosen to test (validate) the accuracy of the AMDM as comparisons can be made with what is already given in the literature.For this, the values of P, K 0 , and nonlinear term in equation ( 16) were set to zero.
e first five natural frequencies, (Ω n(i) , i � 1, . . ., 5) are shown in Table 1 and compared with those obtained by Reference [9].Excellent agreement was found.
e results are shown in Table 2, with, and without, the nonlinear term given in equation ( 16).Here, Ω non(i) are the results with the nonlinear term used and Ω n(i) are the results when the nonlinear term was not used.
It can be seen from Table 2 that the results agree with those of Reference [27] fairly well where the nonlinear term is absent although generally the present results are slightly higher than those calculated by Chen [27].Also shown in Table 2 are the results calculated for large amplitudes. is was effected through the inclusion of the nonlinear term of equation ( 16).
e ratio of results calculated with and without the inclusion of the nonlinear term is also given in Table 2. e presence of the nonlinear terms increases quite substantially the natural frequencies, and it can be seen that this increase grows in line with increasing the mode number.
Figure 3 presents the variation of natural frequency modes with increase in axial force and the foundation stiffness.During the variation of the axial force the elastic stiffness, K 0 , was held constant and for the variation of the elastic stiffness, the axial force, P, was held constant.It is noticeable that the first mode is significantly affected by both variations in comparison with the higher modes, and in particular, the first mode natural frequency is greatly affected by increasing in axial force.
By substituting the converged Ω [n]  n(1) into equations ( 21) and ( 35) and normalizing the result using equation (42), a polynomial can be obtained to describe the first mode shape function.e same procedure can be employed for other natural frequencies to find the mode shapes for higher mode numbers.e variation of the first and third mode shapes for various Winkler parameters is illustrated on Figure 4. Here, the axial force was set to zero.e changes in the mode shape were not too significant until the value of K 0 became reasonably large.It can be seen from the results for the third mode that the increase in elastic stiffness affects both the amplitude and the phase of the shape function.
A similar exercise was carried out to ascertain the trends for the mode shapes when the axial force was increased or when the beam was tensioned as opposed to being compressed.e effects on the first and third modes are shown in Figure 5. Here, the elastic stiffness was held constant at one.It can be seen that the greatest effect was on the first mode and, except at the extrema, increasing the axial load had very little effect on higher modes.
Important to this study was the speed of obtaining accurate converged solutions.An example of the rate of convergence is given in Figure 6 for the case of clamped-free boundary conditions with P � 0 and K 0 � 1.As can be seen, the AMDM method was fast converging with results for the first mode obtained after a very few iterations.Convergence for higher modes took longer, although not prohibitively.
6.2.Clamped-Clamped Uniform Beam.For this case, the boundary conditions are as given in equations ( 17) and (19) with the spring constants becoming κ L0 ⟶ ∞, κ R0 ⟶ ∞, κ L1 ⟶ ∞, and κ R1 ⟶ ∞: Calculations were made for the clamped-clamped beam of the natural frequencies with no axial force and for two values of foundation stiffness as shown in Table 3. e results for the present method compare well with those reported in 6 Shock and Vibration   Shock and Vibration the literature.ese results were obtained with the nonlinear term of equation ( 16) not used.
Table 4 shows results for the first four modes with and without the nonlinear term present.Again Ω non(i) is the result when the nonlinear term is used, and Ω n(i) is the result when the nonlinear term is not used.It can be seen that the present calculated values of Ω n(i) are in good agreement although again most of the results are slightly higher than those calculated by Chen [27].When the nonlinear term is included, the results for the natural frequencies increase with the ratio of increase found to be similar in magnitude to those reported by Mei [30].
Again polynomials were obtained to describe the first three mode shape functions for the clamped-clamped beam.
e first three mode shapes for the clamped-clamped beam are shown in Figure 7 after normalization and compared well with those reported in the literature.
Finally, the present method of calculation was validated using a backbone curve for the clamped-clamped beam as shown in Figure 8. e first mode results reported by Gupta  Table 3: Natural frequencies of the clamped-clamped beam on an elastic foundation with P � 0 and without the nonlinear term.
Reference [28] Reference [29] Present Reference [28] Reference [29] Ω n( 1 [31] were chosen with reasonable agreement found with a similar trend found between the two result sets, as illustrated in Figure 8.Here, the value of the foundation stiffness was set at zero.

Conclusions
A fast, efficient, and accurate method of solution, namely, the Adomian modified decomposition method (AMDM) was developed to calculate natural vibrations of an Euler-Bernoulli beam with large amplitude resting on a Winkler foundation.e method is free of linearization, discretization, perturbation, or a priori assumptions and nonlinear terms are relatively easily treated.A practical advantage of the AMDM is the ease of applying the boundary conditions where the vibrational analysis for different boundary conditions simply involves changing the values of the corresponding parameters with no need to change the solution procedures or the algorithms employed.e numerical comparisons for both boundary conditions used here indicate that the current numerical results are in satisfactory agreement with those found by other methods, with perhaps some of the results found here better than some reported.

Figure 1 :
Figure 1: Beam under axial loading while resting on an elastic foundation.

Figure 2 :
Figure 2: Boundary condition described by rotational and translational flexible ends.

Figure 3 :
Figure 3: Variation of frequency modes with axial force and foundation stiffness.

Figure 4 :
Figure 4: (a) First and (b) third mode shapes for various stiffness parameters with P � 0.

Figure 7 :
Figure 7: First three mode shapes for a clamped-clamped beam with K 0 � 1 and P � 1.

Figure 8 :
Figure 8: Comparison of backbone curves for the 1 st mode of the clamped-clamped beam.

Table 2 :
Natural frequencies of the cantilever beam on an elastic foundation with P � 0 and K 0 � 1.

Table 4 :
Natural frequencies of the clamped-clamped beam on an elastic foundation with P � 0 and K 0 � 1.