On the steady state response of a cantilever beam partially immersed in a fluid and carrying an intermediate mass

This paper presents a study on the nonlinear steady state response of a slender beam partially immersed in a fluid and carrying an intermediate mass. The model is developed based on the large deformation theory with the constraint of inextensible beam, which is valid for most engineering structures. The Lagrangian dynamics in conjunction with the assumed mode method is utilized in deriving the non-linear unimodal temporal equation of motion. The distributed and concentrated sinusoidal loads are accounted for in a consistent manner using the assumed mode method. The non-linear equation of motion is, analytically, solved using the single term harmonic balance (SHB) and the two terms harmonic balance (2HB) methods. The stability of the system, under various loading conditions, is investigated. The results are presented, discussed and some conclusions on the partially immersed beam nonlinear dynamics are extracted.


Introduction
Offshore structures such as piles, oil plat-forms supports, oil-loading terminals and towers surrounded by water is usually modeled as a beam or a column when studying its static or dynamic behavior.Since these structures are relatively flexible due to their high aspect ratio and because they are usually subjected to various excitation loads such as wind loads and wave loads, the prediction of their steady state responses, under various combination of parameters, is extremely needed for design and analysis purposes.
Westergard [1] investigated the hydrodynamic pressure on a rigid dam under earthquake excitation.His investigation resulted in the finding that the magnitude of the hydrodynamic pressure is dependent on the frequency of excitation.In the same direction, many investigations on the hydrodynamic pressure on rigid structures during earthquake excitation can be found in the open literature.Among these investigations are Chopra [2], Chwang [3], Lin [4] and others.When the issue is the natural frequencies of immersed structures, Change and Liu [5] have referred to the findings of some investigations that utilized the transfer matrix approach for approximating the wet-beam deflection.They have used the same approach for a tapered immersed beam that is carrying a tip mass and supported by torsional and linear springs.Their model adopted the Euler-Bernoulli beam theory, with the assumption of small deformations, and the Morsion's condition to account for the fluid effect.Recently, Xing et al. [6] reported the results of their investigation on the natural frequencies of beam water interaction system.They developed a coupled fluid-beam dynamic model in which the fluid domain was modeled by a pressure differential equation and the structure was modeled using twosegments Euler-Bernoulli beam differential equations, with the small deformations assumption.The surface continuity conditions in conjunction with the system fluid and structure boundary conditions were utilized in developing the coupled system solutions.Their calculations showed that the effect of fluid can be taken as an added inertia to the structure.Furthermore, they investigated the effects of free surface wave disturbance and the effects of fluid compressibility.
In the previously cited investigations where the immersed beam dynamics was the concern, only the natural frequencies of a structure that undergoes small deformation were calculated and no attention was given to system nonlinearities that develop as a result of large deformations due to beam slenderness.These nonlin-Shock and Vibration 7 (2000) [179][180][181][182][183][184][185][186][187][188][189][190][191][192][193][194] ISSN 1070-9622 / $8.00  2000, IOS Press.All rights reserved earities are expected to show up in offshore structures due to their inherit flexibility and thus produces amplitude frequency dependence.Recently, Al-Bedoor et al. [7] studied the non-linear natural frequencies of a slender beam partially immersed in water.Their results have documented a parametric study on the frequencyamplitude nonlinear relation under different parameters combinations.The immersed beam system, has shown a hardening amplitude-frequency relation for the first mode and a softening behavior for the second and higher modes.Due to the facts that such offshore structures are extremely flexible and they are usually exposed to different distributed and concentrated loads with harmonically varying characteristics, their steady state response should be studied.
This work is devoted towards studying the nonlinear steady state response of a slender beam partially immersed in water and carrying an intermediate mass.
The model is developed based on the large deformation theory with the constraint of inextensible beam, which is valid for most engineering structures.The Largrangian dynamics in conjunction with the assumed mode method is utilized in deriving the unimodal temporal equation of motion.Moreover, the distributed and concentrated sinusoidal loads are discretized using the assumed mode method.The nonlinear equation of motion are solved using the single term (SHB) and two terms (2HB) harmonic balance method, The selection of this method to analyze the nonlinear equation in the present work is based on the analysis results presented on a wider class of problems in reference [8] which showed that the 2HB can capture the correct behavior of the steady state periodic response of such a system.In addition it was also whown that the 2HB, in comparison with the method multiple scale (MMS), does not generate spurious solutions nor tends to distorte the steady state response curves of such oscillator.Furthermore, in the presnt work the first order stability are presented and discussed.

System description and assumptions
A schematic of the beam under study is shown in Fig. 1.The beam is considered to be uniform of constant length l, cross sectional area A, flexural rigidity EI and density ρ.The beam is vertically mounted, clamped at the base and partially immersed in a fluid up to depth l 1 and carries a lumped mass M at an ar-bitrary intermediate position d along the beam span.The fluid is assumed to be non-viscous, incompressible with a constant density ρ .The thickness of the beam is assumed to be small compared to its length, so the effect of rotary inertia and shear deformation can be ignored.It is assumed also that the beam can only undergo planar felxural vibrations and the amplitude of this vibration may reach any large value, but remains below the limiting value of the slope θ of the elastica |θ| = 90 • .These assumptions are the same that used in studying the planar nonlinear vibrations of a cantilever beam systems [9,10].

Equation of motion
Upon using the coordinate system shown in Fig. 1 for the immersed beam, the elastic potential Energy V e of the beam due to bending is given by; where ζ = s/l is the dimensionless arc length and R(ζ) is the of curvature of the beam neutral axis.The exact of curvature R(ζ) takes the form [9], where λ = 1/l and the prime denotes differentiation with respect to the non-dimensional arc length ζ.The system potential energy must contain the gravitational potential energy develops as a result of axial shortening due to transverse deformations.An expression of the axial shortening based on the assumption that no axial deflection [11], is given by; The gravitational potential energy of the system, upon using equation (3), can be written as; where M is the lumped mass located at the position η = d/l.The resulted potential energy of the system can obtained by adding equation ( 4) to (1), i.e.; The potential energy of the system can be expressed in terms of the beam transverse deflection variable y, by noting that the variables x and y are related through the subsidiary relation [9] Then by using equation (6) and its derivative to eliminate x and x from equation (5), upon expanding the term (1 − λ 2 y 2 ) 1/2 into a power series and noting that (λ 2 y 2 ) < l, equation (5) becomes when the nonlinear terms retained up to fourth order; where µ = M/ρAl.Next the kinetic energy T of the system is given by; where C m is the inertia coefficient of the additional mass of the fluid and K m = ρ /ρ as imposed by Chang and Liu [5].In order to express the system kinetic energy expressions in terms of the beam transverse deflection variable y, it is to be noted also, that the beam under study is assumed to be inextensible, which implies that the length of the neutral axis of the beam remains constant during vibration.This imposes the relation [12] (1 + λx ) 2 + (λy After mathematical manipulations and simplifications, i.e. expanding into power series, noting that (λ 2 y 2 ) < 1, and integrating from 0 to an arbitrary value of ζ, equation ( 9) can be rewritten in the form; Differentiating with respect to time t, squaring and retaining the non-linear terms up to fourth order leads to Upon substituting equation (11) into equation ( 8), system kinetic energy T of the considered beam system becomes; where ζ 1 = l 1 /l, using equations ( 7) and ( 12) one obtains, after factoring out ρAl/2, the system Lagrangian L, L = T − V , as where β 2 = EIλ 3 /ρAl is the linear frequency parameter and g 0 = ρAl 3 g/EI is the dimensionless gravity parameter.
To obtain an approximate ordinary differential equation in time t, the assumed mode method is utilized, such that, where q(t) is an unknown time modulation of the assumed mode shape and φ(ζ) is the normalized mode shape function of the linear cantilever beam, i.e. l 0 φ 2 (ζ)dζ = 1, which is assumed to remain selfsimilar (i.e.independent of motion amplitude) during the motion.In this work Galerkin's method is used, whereby φ(ζ) is the eigenfunction of the n − th mode of the cantilever beam, given in many vibration text books, as Where p is the n− th roots of the frequency equation 1 + cos p cos hp = 0, the first four roots of equation are, 1.875104, 4.694091, 7.854757 and 10.99554.
Substituting equation ( 14) into equation ( 13), one obtains the discrete unimodal beam Lagrangian, which can be expressed as where α i (i = 1, 4) are dimensionless coeffeicints associated with the selected mode shapes used to discretize the system Lagrangian, and are defined as follows: To study the forced planar response of the beam system, two types of periodic excitations loads are used the first is a concentrated load F 0 acting at an arbitrary ζ c along the span of the beam and the second is a distributed load, given as F 0 (s) cos(Ωt).Both loads are assumed to act only in the y direction, i.e. beam transverse direction.Upon the application of the Euler-Lagrange equation, where Q is the generalized force, which can be determined from the principle of virtual work method, δW = Q • δy, one obtains the discrete beam nonlinear equation as follows, where α 5 = φ(ζ c ) for the concentrated load and α 5 = 1 ζ1 φ(χ)dχ for the distributed load.It is to noted that some of the coefficients α i in equations (17-20), increase sharply and attain relatively large values at the higher modes of the beam.Therefore, for convenience, equation ( 22) is scaled and converted to the dimensionless form Where dots are now derivatives with respect to the dimen- l is the dimensionless displacement amplitude at the point of maximum deflection, and p 2 = ωβ is the dimensionless frequency, and ω is the frequency of the assumed mode of the associated linear cantilever beam.
For some stability analysis, structure and water effective damping is assumed to be viscous, with damping coefficient δ, which can be added to the equation of motion to take the form; Equation ( 25) belongs to the same class of nonlinear oscillators studied by Al-Qaisia and Hamdan [8].This equation describes the nonlinear planar flexural motion of the in-extensible beam shown in Fig. 1.In this equation the terms ε 1 u 2 ü and ε 1 u u2 are inertial nonlinearities due to kinetic energy of axial motion and arise in this equation as a result of using the in-extensibility condition and its derivatives given in equations ( 10) and (11).The first of these nonlinear terms has a softening effect (i.e.leads to a decreasing natural frequency with decreasing motion ampitude, or equivalently tends to bend the frequency response curve to the left), while the second nonlinear term has a hardeneing effect (i.e.leads to an increasing natural frequency with increasing vibration amplitude, or equivalently tends to bend the frequency response curve to the right).The last of the nonlinear term ε 2 u 3 , in equation ( 25), is a hardening static type due to potential energy stored in bending and arise in this equation as a result of using nonlinear curvature.Thus, the response of the above nonlinear oscillator is controlled by two competing sofetning and hardening nonlinearities which exhibits fundamnetaly two different response characterestics [8,13].As to which of these nonlinearities dominate the response will depend on the relative value of the coefficints ε 1 and ε 2 in this equation, which are functions of the sysetm parameters and the assumed mode shape.Based on the results presented on [10,13], the steady state response of this oscillator is expected to be of the softening type, when ε 1 /ε 2 1.6, otherwise of the softening type.Generally speaking, in the present problem and for the first mode ε 1 /ε 2 < 1.6, and for the second and higher modes ε 1 /ε 2 1.6 and increases as the mode number increase.Thus one expects that the beam response at the first mode to be of the hardening type while for the second and higher modes this response is of softening type.Approximations to the steady state for the first mode and for the second mode are presented and discussed in the next section.

Solutions of the nonlinear model
Approximate analytical solutions for the periodic steady state response, having the same period as the excitation, of the nonlinear, single mode temporal equation of motion (25) of the beam system described in Section 2, are presented.These solutions are obtained using single term harmonic balance method (SHB) and two terms harmonic balance method (2HB).A new time T = Ωτ is first introduced so that equation (25) becomes; where dots are T derivative and the unknown phase φ has been added to the excitation so that one can obtain a fundamental harmonic response containing only one trigonometric term.

Single term harmonic balance (SHB)
According to the SHB method, an approximate solution of equation ( 26), takes the form; where A is the steady state response amplitude.Substituting equation (27) into equation (26), neglecting third harmonics which arise, and equating coefficients of first harmonics, one obtains the following equations: where P = ε 3 F , the steady state frequency response is obtained by squaring and adding equations ( 28) and (29) and solving for Ω 2 as a function of A; this yields where Second Mode, P = 0.2, δ = 0.02 Equation (30), yields two real solutions for Ω provided that the radical term is real and less than R 1 ; a single real solution is obtained when the radical term is zero or greater the R 1 , and no real solution exists when R 2  1 − R 2 < 0. The steady state frequency response curves obtained using equation (30) are presented and discussed in next section.

Two terms harmonic balance (2HB)
An improved harmonic balance solution can be obtained by including higher harmonics in the assumed solution in equation (27).It is to be noted that the problem of selecting the 'right' combination of the leading harmonics form Fourier series approximation to the response which will lead to the "correct" behavior of the predicted response becomes a difficult task, especially when the nonlinearities are not relatively small.In the present work, the two terms harmonic balance (2HB) solution is selected based on a previous finding by Al-Qaisia and Hamdan [8].It was shown that the 2HB, can predict the correct behavior of the steady state frequency response and an appreciable improvement of the accuracy of the predicted responses can be obtained even when the nonlinearity is relatively strong.In this work, only one more term is added to equation ( 27), whereby the two-terms approximation,having the same period as the excitation, to the steady state solution of the system in equation ( 26) with odd nonlinearities takes the form; Substituting equation (33) into equation (26) and using the same procedure followed previously in Section 3.1 and neglecting the higher-order harmonics, one obtains the following coupled nonlinear algebraic equations for A 1 , A 3 , B 3 and the phase φ; These equations may be expressed in a more convenient form as follows.First, squaring and adding equations (34) and (35) and solving for Ω, leads to Where: Next, equations ( 36) and (37) are solved implicitly for A 3 and B 3 , respectively: Equation (38) can be written using the form where R 3 and R 4 , can be calculated from equation (39), so that, R 3 = (−b/2a) and R 3 = (c/a).Equation ( 42) has two real solutions provided that , and no real solution exists when R 2  3 < R 4 .Equations ( 40), ( 41) and (42) were solved iteratively with an accuracy of 10 −6 to define steady state solution.The steady state frequency response curves obtained using these equations are, for convenience, presented and discussed in Section (4).

Stability of the steady state response
The stability of the steady state response of the fundamental harmonic approximation (27), is examined by introducing a small perturbation ν(T ), i.e. by substituting u(T ) = A cos(T )+ν(T ), into the equation ( 26), followed by use of the steady state conditions (28) and (29).This leads to the following non-linear variational equation; The stability is governed by the linear version of equation ( 43).In addition, the excitation term on the right hand side of equation ( 43) is deleted, because it has no influence on the stability of the response ν(T ).The linear stability is governed by the standard form of the damped Mathieu equation. Where: Second Mode, P=0.2 , δ = 0.02 The approximation of the stability boundaries of equation ( 44) associated with the principal parametric resonance is given by [14].45) are such that the point α * (q * ) lies between the curves (46), the steady state solution is unstable to small disturbances.Thus the conditions for instability may be stated in term of A, Ω, δ, ε 1 , ε 2 as follows: The results obtained for the stability of the steady state response of the fundamental harmonic approximation (27), using the instability conditions (47) and (48), are presented in the next section.

Results and discussion
The steady state frequency response of the nonlinear, single mode, temporal equation of motion (26) of the immersed cantilever beam carrying an intermediate lumped mass "Fig.1", was calculated, for various values of system parameters ε 1 , ε 2 , δ, ε 3 and excitation level P , analytically by using the single term and two terms harmonic balance methods (SHB and Table 1 Values of the parameters in the temporal equation (25) of the beam system shown in Fig. 1 with physical characteristics as described in Section (4), for the first four linear eigenfunctions of the cantilever beam.Where the subscript c corresponding to the concentrated load and the subscript d corresponding to the distributed load.and length l = 2 m.The fluid is taken to be water ρ = 1000 Kg/m 3 , so that K m = 1/2.8.The added inertia coefficient was taken as in [5], C m = 1.0, and the non-dimensional gravity parameter g 0 is calculated to be 0.015.To calculate the parameters ε 1 , ε 2 and ε 3 in the single mode nonlinear temporal equation ( 25 As mentioned in the previous section P = ε 3 F , where ε 3 is equal to ε 3c of the first mode of vibration.Consequently, the value of F = F 0 /ρAl 2 β 2 , is calculated for the distributed sinusoidal load and for the higher modes.The calculated values of α i and the parameters ε 1 , ε 2 and ε 3 for the two types of excitation loads and the immersed and dry beam system are presented, for convenience, in Table 1. Figure 2 shows the steady state frequency response of the immersed beam governed by equation (26) using single term harmonic balance (SHB) and two terms harmonic balance (2HB) methods for the first mode and for both cases, the concentrated load and the distributed load.These results show that for the selected system parameters, the first mode steady state frequency response exhibits, as expected, a hardening behavior, due to the fact that (ε 1 /ε 2 ∼ = 1.12) < 1.6 [8,13].The steady state frequency response of the second mode, was obtained also using (SHB) and (2HB) methods for the two types of excitation loads, is shown in Fig. 3.The steady state frequency response in this case exhibits a softening behavior as the ratio of (ε 1 /ε 2 ∼ = 4.6) > 1.6, i.e. the response is dominated by the inertia nonlinearities.
These results show that the steady state response curve bend to the right, i.e. exhibits a hardening behavior, but with an increasing slope as the driving frequency is increased from a value below the lin- ear resonance frequency.At relatively moderate values of response amplitude the slope of the response curve, for the lightly damped system, begins to increase rapidly whereby the response curve exhibits a relatively "sharp" resonance peak similar to that of a linear oscillator.This behavior of the systems with competing inertia (softening) and static (hardening) nonlinearities is due to the fact that as the response amplitude becomes relatively large (for the lightly damped system) the backbone curve (i.e. the curve defining the free response amplitude-frequency variation) becomes nearly constant independent of the response amplitude, i.e. similar to that of a linear oscillator [13].
The results in Fig. 2 also show that the SHB solution in comparison with the 2HB solution may produce a significant "shift" to the right in the steady state frequency response curves.These results as well as the more detailed ones in [8,13] indicate that the SHB solution or equivalently the first order perturbation solu-tions [8] may not display the correct frequency response behavior of these types of oscillators.
It is to be noted also, from Fig. 3, that the steady state frequency responses, of the second mode for the distributed excitation load, obtained by SHB and 2HB are of the same order and the amplitude is small in comparison that of the concentrated excitation load.This is due to the fact that the value of ε 3c /ε 3d ∼ = 15.4,which means that the excitation level is very small.Figures 4 and 5 show the first order instability boundaries obtained using equations ( 47) and (48) and the steady state frequency response curves obtained using SHB for the first and second modes with are shown as well as the frequency response curves obtained by using SHB for the first mode with p = 0.2 and δ = 0.02.An important feature of these results is that for both the hardening (first mode, Fig. 4) and the softening (second mode, Fig. 5) these instability curve do not intersect solution steady state response curves "SHB" at the ver-Second Mode, P = 0.2, δ = 0.02 tical tangency point as one may expect.As mentioned in Subsection (3.3), the instability conditions (47) and ( 48) is governed by A, Ω, δ, ε 1 , ε 2 and to show the effect of the ratio ε 1 /ε 2 on the instability boundary, decreasing the ratio of ε 1 /ε 2 form 1.12 to 0.85 for the same oscillator, one can see form Fig. 6 that the instability curve intersects the steady state response curve at the point with vertical tangency.This result indicate that, there is a need for a farther more in depth stability analysis for better understanding of the instability mechanism in this type of oscillators.
Figure 7 shows the steady state frequency responses for the first mode obtained by using the SHB, for different excitation levels P = 0.2, 0.5 and 1.0.These results show that, as the excitation increases, at higher amplitudes, the system reveals a linear behavior as the excitation level increases.
In Figs 8 and 9, the steady state frequency responses using SHB are shown for the first and second mode, respectively, for the wet beam case and the dry beam case.As one, can expect, the effect of the additional mass of fluid (K m = ρ /ρ) is clear on the responses, i.e., the nonlinear natural frequency is decreased for the immersed beam compared to the dry beam, as one can see from the backbone curves.

Conlusions
The present work studied the steady state frequency response of a slender cantilever beam partially immersed in water and carrying an intermediate mass.
The assumption of the inextensibility condition were taken as a means to account the inertia nonlinearities, which have a strong influence on the steady state frequency response curves of the beam system.The nonlinear equation of motion was derived using Euler-Lagrange method in conjunction with the assumed mode method.The steady state responses under the effect of sinusoidal distributed and concentrated loads were obtained.The steady state frequency response curves are presented for the first two modes, which are obtained using the method of harmonic balance, SHB and 2HB.Appreciable improvement of the accuracy of the predicted response was obtained using two terms harmonic balance 2HB.It appear that the results presented in this work indicate that first order approximate solutions and stability results may not adequately uncover the actual behavior of this type of oscillator even when the system is relatively weakly nonlinear, i.e. ε 1 and ε 2 are small compare to unity.The steady state responses show that the effect water "added inertia" is to decrease the natural frequency as related to the vibration amplitude as shown by the left skewed backbone curves.

Aknowlegment
The authors would like to acknowledge the supports of both the University of Jordan, Amman-Jordan and King Fahd University of Petroleum and Minerals, Dhaharan-Saudi Arabia.

Fig. 1 .
Fig. 1.A schematic of the immersed cantilever beam under consideration.

P=0. 2 ,Fig. 6 .Fig. 7 .
Fig. 6.Same as Fig. 4, but for ε 1 = 0.85 and ε 2 = 1.0.2HB),equations (30) and (38), respectively.As an example, the beam if Fig.1is taken to be an Aluminum ), for selected values of the attached mass magnitude ratio µ and relative position η, relative position of the applied concentrated load ζ c and relative depth ζ 1 of the wet section, the integrals in equation (17-20) defining the coefficients were evaluated numerically by using the symbolic manipulator program Derive.The selected parameters values used in this work were µ = 0.25, η = 0.75, ζ c = 0.7, ζ 1 = 0.5, δ = 0.02 and the excitation level P is selected to be within the range 0.2-1.0.