- Email: [email protected]
Dynamic response of continuous thick beams subjected to a moving force
Journal of'Sound and Vibration (1986) 109(1), 169-173
LETTERS T O THE EDITOR DYNAMIC RESPONSE OF CONTINUOUS THICK BEAMS SUBJECTED TO A MOVING FORCE 1. INTRODUCTION The solution of continuous thick beams subjected to a moving load is considered through the use of Timoshenko beam theory which takes into account the effect of transverse shear deformation and rotary inertia. The purpose of this note is to describe a numerical approach for solving this problem. In this approach, the solution is formulated as a series of products of the beam eigenfunctions, determined by using the "dynamic flexuralstiffness matrix approach", and time dependent functions as detailed below. 2. NUMERICAL PROCEDURE
2.1. Timoshenko's beam equations A basic assumption of Timoshenko beam theory is that the cross-sectional rotation, ~, depends not only on the slope of the (total) deflection curve, as in Bernouilli-Euler theory, but also on the value of the shear strain. Thus the deformation of any point in the beam is related to the value__ssof two reference quantities, w (the deflection) and ~b, at the neutral axis (rather than w alone as as in Bernouilli-Euler theory). It follows that the dynamic Timoshenko beam problem is governed by a pair of coupled differential equations of the form [1] \OX
--P Ot2 -
'
\OX
-pA-~-~=O.
(1)
Here E is the modulus of elasticity, G the modulus of rigidity, I the second moment of area, ,X, the area of the cross-section, p the material mass density and K 2 the shear factor. On the assumption of harmonic motion, i.e.
w(x, t)= W(x) e iw,,
~b(x, t) = ~ ( x ) e iw',
(3)
the general solution of the governing differential equations for the flexural modes of the vibration of a uniform beam of length l is [2, 3] W(sc) = C1 sin fl~:+ C2 cos fl~+ (?3 sinh a ~ + (?4 cos a~:,
(4)
g'(~) = Clkl cos [3~- C2kl sin fl~+ C3k2 cosh ag+ C4k2 sinh a~:,
(5)
kl_b(fl2-s2) fll '
(6)
where
k2=b(a2+s2) al '
r2 = I /~i~'
s2 -
E1 G.Z,------1" K2
It is assumed here that bZr2s2<~ 1 (b2r2s2> 1 is treated in Appendix I of reference [4]). to is the radian frequency and Ct, C2, C~, and (?4 are constants which may be found from equations (4) and (5) with the following boundary conditions: at~:=0W=Wa
and
gt=9-o,
at~:=l W=Wb
and
g'=9"b.
(7)
169
0022-460x/86/160169+05 $03.00/0
9 1986 AcademicPress Inc. (London) Limited
170
LETTERS TO T H E E D I T O R
0I
l
b iX
I
T~ v b .
1
1 Y
Figure 1. Forces and displacements at the ends of an element.
2,2. Dynamic flexural-stiffness matrix The displacements, cross-sectional rotations, applied forces and moments at the ends a and b of a uniform beam of length I are shown in Figure 1. Thus, the force-displacement relationship can be expressed in matrix.form as follows [4]:
M4=ls2 s5-s,
soll 4
vq
llW l '
(8)
Sl = (EIIl3)(b2/ot/3y)(a + r#3)(sinh a/3 cos/3 + r/sin/3 cosh a), $2 = ( E I / 1 2 ) ( k J / y ) [ ( a - ~7/3)(cosh a cos/3 - 1) + (/3 + r/a) sin/3 sinh a], $3 = ( - E I /13)( b 2/ a/3y )( a + r/fl)(sinh a + r/sin/3), $4 = ( E I / IE)( rl / y)( a 2+ /32)( cosh a - c o s / 3 ) , Ss = ( E I / l ) [ ( c t + r//3)/y](sin/3 cosh ot - r/cos/3 sinh a), $6 = ( E I / l ) [ ( a + r/fl)/y] (rl sinh a - s i n / 3 ) , y = 2r/(1 - c o s h cr cos/3) + (1 - r t 2) sin/3 sinh a,
7/= kl/k2.
(9)
Systems that consist of several beams may be assembled from the dynamic member stiffness matrices of all individual members in the usual manner. Thus, if there are n independent non-zero displacements at the joints and the boundaries of the system, we obtain n equilibrium equations, which can be written in terms of displacements d as Sd=P,
(I0)
where S is an n x n matrix. The elements of S depend upon the frequency to. For free vibration problems the load vector P is zero and the natural frequencies are determined from det ISl =0.
(11)
The eigenvalues to of equation (11) may be obtained by applying the method of regula falsi [5]. 2.3. Eigenvalues and eigenfunctions Corresponding to each eigenvalue w, there is an eigenvector d defined by equation (10). The elements of the vector d represent all the degrees of freedom: i.e., the rotations
171
LETTERS TO THE E D I T O R
and the displacements at the joints. These can be determined by setting any element of the vector d in equation (10) equal to 1, and finding the relative values of the other element. For a continuous beam the eigenfunctions take the following form of equations (4) and (5): W.. = (Ci)., sin/3~r~+ (C2)n, cos/3.,~+ (C3)., sinh a~r~+ (C4)~r cosh a.~,
~n; = (Clkl)... cos f l . ~ -
(12)
(C2kl)~, sin fl.,lj+ (C3k2),.rcosh a.~:+ (C4k2)., sinh ct.~:.
(13)
Here W., and 9'., are the deflection and the cross-sectional rotation of the rth element for the nth mode; (a).r, (fl).,. (kl),.r and (k2),.. can be determined from equation (6). The constants C~. C2. C3 and (74 can be determined by substituting of the appropriate deflections and the cross-sectional rotations obtained from the eigenvector d of equation (10) into equations (12) and (13) as follows:
Wo]( ~. Wb
~b.,
0
=
1
0
1
\r,,=,,.,.1 ~J(C2)., [
(k.)~,
0
(k2).,
0
sin fin, \(kl)~r COSfin.
cos ft., -(kl)~,sin/3..
sinh a., (k:).r cosh am,
cosh a.r
~
I I (C3)"~1" (k2)nrsinha.4[(C4).r j (14)
2.4.
Energyformulation
To determine the response of the continuous beam shown in Figure 2, the dynamic deflection w(x, t) and the cross-sectional rotation ~(x, t) may be expressed by the summation of the modal components: p
p
w(x, t)= • ff'.(x)" An(t),
~(x, t)= • gr.(x)" Bn(t).
n=l
A
I.
It
(15,16)
n=l
.A.
"~
lz
LI. . . . .
I
'-
'
L
'1 I
Figure 2. A N - s p a n continuous beam subjected to a moving load P.
Here ~/. and ~ . are the deflection and the cross-sectional rotation of the nth mode of the whole structure and An(t) and B.(t) are the time dependent modal amplitudes for the deflection and rotation respectively. The strain energy, U, of the system described by Timoshenko may be expressed as 1
i.
2
U = ~ Io [EI(~x~) z+K
-
aw
GA("~x-lp)
2
]dx,
(17,
where L is the total length of the continuous beam. Thekinetic energy, T, of the vibrating beam is, where (') denotes d( )/dt,
1/[ [p,4(k)2+pI(~t)2]dx.
T=~
(18)
172
LE'I'FERS
TO
THE
EDITOR
By using equations (15) and (16) with equation (17), the strain energy may be re-expressed as
1
U=~CTKC, where
Cr={A,(t)
B,(t)
A2(t)
""
As(t)
(19) /~,(t)
...
Ap(t)
BAt)},
and K is the stiffness matrix. Similarly, by substituting equations (15) and (16) into equation (18) the kinetic energy, T, may be re-written as (20)
T=I crMC.,
where M is the mass matrix. 2.5. Dynamic response to a moving force When a force P moves at constant speed v, on the sth span, as shown in Figure 2, the external work done, We, is
~,
Vts
Applying Lagrange's equation leads to a set of simultaneous differential equations of the second order, Mi~ + KC = Q, where
Qr=P{W~(vtJl~)
0
""
W~s(vt,/l~)
(22) 0
""
Wp,(vt,/l~)
0}.
Note that t~ = 0 is the time when the load reaches the left node of the sth element if the load is moving to the right. Equation (22) can be solved exactly [6] or numerically. One numerical technique used to obtain a solution is an implicit finite difference scheme known as Houbolt's [7, 8]. For the case given below this method is stable irrespective of the value of the time interval At [91. 2.6. Houbolt's method For any point j on the time co-ordinate, t, (0 <. ts<~ l,/u), the first and the second time derivative may be expressed as polynomials of the third degree passing through points j , j - 1 , j - 2 and j - 3 of the grid, as follows [8]: C ) = [1/(At)2](2Cj - 5Cj_~ + 4Cj_2 + Cj-3), C'j = [1~6At](11Cj
--
18Cj_, + 9Cj_2 - 2Cy_3).
(23)
Thus, equation (22) can be replaced by the expression (K+ [2/(flt)2]M)Cj = [1/(At)2]M(5Cj_~ - 4Cj_2 + C~-3) + Qj.
(24)
It is clear that if the solution is known for the three preceeding points in time, the solution at point j can be obtained by solving equation (24). In the first two steps of the procedure not all the data required for evaluating the right hand side on equation (24) are available. However to overcome this shortage, the Newmark method [8], in which constant average accelerations are assumed, was employed and the final algorithm reads as follows.
LETTERS TO THE EDITOR
173
F o r j = 1, 2(t, =jAr) KCj = M{[4/(At)2]Cj_I + [4/At]C~_~ + C"j-a} + Qj C " / = [ 4 / ( A t) 2] (Cj - Cj_I) - [4/At] C's_ 1- C"~_I,
(25)
C'j = C'j_, + CAt~2)(C'j_, + C3). For j >i 3 { K + [2/(At)2]M}Cj = [1/(A t)2]M(5Cj_~ - 4Cj_2 + Cj_3) + Qj.
(26)
The initial vectors Co and C'o for any intermediate element can be obtained from the end conditions o f the preceeding element. Finally, substituting the values of the nodal amplitudes A'~s and B ' s calculated from equations (25) and (26) into equations (15) and (16), the displacements and the rotati6ns o f any point along the b e a m at any time can be obtained.
Department of Civil Engineering, New Mexico State University, Las Cruces, New Mexico 88003, U.S.A.
O. L. ROUFAEIL
(Received 31 July 1985, and in revisedform 14 February 1986) REFERENCES
1. S. TIMOSHENKO 1921 Philosophical Magazine 41, 288-289. On the correction for shear of the differential equation for transverse vibration of prismatic bars. 2. R. W. TRAILL-NASH and A. R. COLLAR 1953 Quarterly Journal of Mechanics and Applied Mathematics 6, 186-222. The effect of shear flexibility and rotary inertia on bending vibration of beams. 3. T. C. HUANG 1961 Journal of Applied Mechanics 28, 579-584. The effect of rotary inertia and shear .deformation on the frequency of normal mode equations of beam with simple end conditions. 4. W. P. HOWSON and F. W. WILLIAMS 1973 Journal of Sound and Vibration 26, 503-515. Natural frequencies of frames with axially loaded Timoshenko members. 5. B. WENDORFF 1966 Theoretical Numerical Analysis. New York. Academic Press. 6. G. B. WARBURTON 1976 The Dynamical Behaviour of Structures. Oxford: Pergamon Press, second edition. 7. J. C. HOUBOLT 1950 Journal of Aeronautical Science 17, 540-550. A recurrence matrix solution for the dynamic response of elastic air craft. 8. K. J. BATHE 1982 Finite Element Procedure in Engineering Analysis. Englewood Cliffs, New Jersey: Prentice-Hall. 9. D. E. JOHNSON 1966 American Institute of Aeronautics and Astronautics Journal 4, 1450-1451. A proof of the stability of the Houbolt method.
LETTERS T O THE EDITOR DYNAMIC RESPONSE OF CONTINUOUS THICK BEAMS SUBJECTED TO A MOVING FORCE 1. INTRODUCTION The solution of continuous thick beams subjected to a moving load is considered through the use of Timoshenko beam theory which takes into account the effect of transverse shear deformation and rotary inertia. The purpose of this note is to describe a numerical approach for solving this problem. In this approach, the solution is formulated as a series of products of the beam eigenfunctions, determined by using the "dynamic flexuralstiffness matrix approach", and time dependent functions as detailed below. 2. NUMERICAL PROCEDURE
2.1. Timoshenko's beam equations A basic assumption of Timoshenko beam theory is that the cross-sectional rotation, ~, depends not only on the slope of the (total) deflection curve, as in Bernouilli-Euler theory, but also on the value of the shear strain. Thus the deformation of any point in the beam is related to the value__ssof two reference quantities, w (the deflection) and ~b, at the neutral axis (rather than w alone as as in Bernouilli-Euler theory). It follows that the dynamic Timoshenko beam problem is governed by a pair of coupled differential equations of the form [1] \OX
--P Ot2 -
'
\OX
-pA-~-~=O.
(1)
Here E is the modulus of elasticity, G the modulus of rigidity, I the second moment of area, ,X, the area of the cross-section, p the material mass density and K 2 the shear factor. On the assumption of harmonic motion, i.e.
w(x, t)= W(x) e iw,,
~b(x, t) = ~ ( x ) e iw',
(3)
the general solution of the governing differential equations for the flexural modes of the vibration of a uniform beam of length l is [2, 3] W(sc) = C1 sin fl~:+ C2 cos fl~+ (?3 sinh a ~ + (?4 cos a~:,
(4)
g'(~) = Clkl cos [3~- C2kl sin fl~+ C3k2 cosh ag+ C4k2 sinh a~:,
(5)
kl_b(fl2-s2) fll '
(6)
where
k2=b(a2+s2) al '
r2 = I /~i~'
s2 -
E1 G.Z,------1" K2
It is assumed here that bZr2s2<~ 1 (b2r2s2> 1 is treated in Appendix I of reference [4]). to is the radian frequency and Ct, C2, C~, and (?4 are constants which may be found from equations (4) and (5) with the following boundary conditions: at~:=0W=Wa
and
gt=9-o,
at~:=l W=Wb
and
g'=9"b.
(7)
169
0022-460x/86/160169+05 $03.00/0
9 1986 AcademicPress Inc. (London) Limited
170
LETTERS TO T H E E D I T O R
0I
l
b iX
I
T~ v b .
1
1 Y
Figure 1. Forces and displacements at the ends of an element.
2,2. Dynamic flexural-stiffness matrix The displacements, cross-sectional rotations, applied forces and moments at the ends a and b of a uniform beam of length I are shown in Figure 1. Thus, the force-displacement relationship can be expressed in matrix.form as follows [4]:
M4=ls2 s5-s,
soll 4
vq
llW l '
(8)
Sl = (EIIl3)(b2/ot/3y)(a + r#3)(sinh a/3 cos/3 + r/sin/3 cosh a), $2 = ( E I / 1 2 ) ( k J / y ) [ ( a - ~7/3)(cosh a cos/3 - 1) + (/3 + r/a) sin/3 sinh a], $3 = ( - E I /13)( b 2/ a/3y )( a + r/fl)(sinh a + r/sin/3), $4 = ( E I / IE)( rl / y)( a 2+ /32)( cosh a - c o s / 3 ) , Ss = ( E I / l ) [ ( c t + r//3)/y](sin/3 cosh ot - r/cos/3 sinh a), $6 = ( E I / l ) [ ( a + r/fl)/y] (rl sinh a - s i n / 3 ) , y = 2r/(1 - c o s h cr cos/3) + (1 - r t 2) sin/3 sinh a,
7/= kl/k2.
(9)
Systems that consist of several beams may be assembled from the dynamic member stiffness matrices of all individual members in the usual manner. Thus, if there are n independent non-zero displacements at the joints and the boundaries of the system, we obtain n equilibrium equations, which can be written in terms of displacements d as Sd=P,
(I0)
where S is an n x n matrix. The elements of S depend upon the frequency to. For free vibration problems the load vector P is zero and the natural frequencies are determined from det ISl =0.
(11)
The eigenvalues to of equation (11) may be obtained by applying the method of regula falsi [5]. 2.3. Eigenvalues and eigenfunctions Corresponding to each eigenvalue w, there is an eigenvector d defined by equation (10). The elements of the vector d represent all the degrees of freedom: i.e., the rotations
171
LETTERS TO THE E D I T O R
and the displacements at the joints. These can be determined by setting any element of the vector d in equation (10) equal to 1, and finding the relative values of the other element. For a continuous beam the eigenfunctions take the following form of equations (4) and (5): W.. = (Ci)., sin/3~r~+ (C2)n, cos/3.,~+ (C3)., sinh a~r~+ (C4)~r cosh a.~,
~n; = (Clkl)... cos f l . ~ -
(12)
(C2kl)~, sin fl.,lj+ (C3k2),.rcosh a.~:+ (C4k2)., sinh ct.~:.
(13)
Here W., and 9'., are the deflection and the cross-sectional rotation of the rth element for the nth mode; (a).r, (fl).,. (kl),.r and (k2),.. can be determined from equation (6). The constants C~. C2. C3 and (74 can be determined by substituting of the appropriate deflections and the cross-sectional rotations obtained from the eigenvector d of equation (10) into equations (12) and (13) as follows:
Wo]( ~. Wb
~b.,
0
=
1
0
1
\r,,=,,.,.1 ~J(C2)., [
(k.)~,
0
(k2).,
0
sin fin, \(kl)~r COSfin.
cos ft., -(kl)~,sin/3..
sinh a., (k:).r cosh am,
cosh a.r
~
I I (C3)"~1" (k2)nrsinha.4[(C4).r j (14)
2.4.
Energyformulation
To determine the response of the continuous beam shown in Figure 2, the dynamic deflection w(x, t) and the cross-sectional rotation ~(x, t) may be expressed by the summation of the modal components: p
p
w(x, t)= • ff'.(x)" An(t),
~(x, t)= • gr.(x)" Bn(t).
n=l
A
I.
It
(15,16)
n=l
.A.
"~
lz
LI. . . . .
I
'-
'
L
'1 I
Figure 2. A N - s p a n continuous beam subjected to a moving load P.
Here ~/. and ~ . are the deflection and the cross-sectional rotation of the nth mode of the whole structure and An(t) and B.(t) are the time dependent modal amplitudes for the deflection and rotation respectively. The strain energy, U, of the system described by Timoshenko may be expressed as 1
i.
2
U = ~ Io [EI(~x~) z+K
-
aw
GA("~x-lp)
2
]dx,
(17,
where L is the total length of the continuous beam. Thekinetic energy, T, of the vibrating beam is, where (') denotes d( )/dt,
1/[ [p,4(k)2+pI(~t)2]dx.
T=~
(18)
172
LE'I'FERS
TO
THE
EDITOR
By using equations (15) and (16) with equation (17), the strain energy may be re-expressed as
1
U=~CTKC, where
Cr={A,(t)
B,(t)
A2(t)
""
As(t)
(19) /~,(t)
...
Ap(t)
BAt)},
and K is the stiffness matrix. Similarly, by substituting equations (15) and (16) into equation (18) the kinetic energy, T, may be re-written as (20)
T=I crMC.,
where M is the mass matrix. 2.5. Dynamic response to a moving force When a force P moves at constant speed v, on the sth span, as shown in Figure 2, the external work done, We, is
~,
Vts
Applying Lagrange's equation leads to a set of simultaneous differential equations of the second order, Mi~ + KC = Q, where
Qr=P{W~(vtJl~)
0
""
W~s(vt,/l~)
(22) 0
""
Wp,(vt,/l~)
0}.
Note that t~ = 0 is the time when the load reaches the left node of the sth element if the load is moving to the right. Equation (22) can be solved exactly [6] or numerically. One numerical technique used to obtain a solution is an implicit finite difference scheme known as Houbolt's [7, 8]. For the case given below this method is stable irrespective of the value of the time interval At [91. 2.6. Houbolt's method For any point j on the time co-ordinate, t, (0 <. ts<~ l,/u), the first and the second time derivative may be expressed as polynomials of the third degree passing through points j , j - 1 , j - 2 and j - 3 of the grid, as follows [8]: C ) = [1/(At)2](2Cj - 5Cj_~ + 4Cj_2 + Cj-3), C'j = [1~6At](11Cj
--
18Cj_, + 9Cj_2 - 2Cy_3).
(23)
Thus, equation (22) can be replaced by the expression (K+ [2/(flt)2]M)Cj = [1/(At)2]M(5Cj_~ - 4Cj_2 + C~-3) + Qj.
(24)
It is clear that if the solution is known for the three preceeding points in time, the solution at point j can be obtained by solving equation (24). In the first two steps of the procedure not all the data required for evaluating the right hand side on equation (24) are available. However to overcome this shortage, the Newmark method [8], in which constant average accelerations are assumed, was employed and the final algorithm reads as follows.
LETTERS TO THE EDITOR
173
F o r j = 1, 2(t, =jAr) KCj = M{[4/(At)2]Cj_I + [4/At]C~_~ + C"j-a} + Qj C " / = [ 4 / ( A t) 2] (Cj - Cj_I) - [4/At] C's_ 1- C"~_I,
(25)
C'j = C'j_, + CAt~2)(C'j_, + C3). For j >i 3 { K + [2/(At)2]M}Cj = [1/(A t)2]M(5Cj_~ - 4Cj_2 + Cj_3) + Qj.
(26)
The initial vectors Co and C'o for any intermediate element can be obtained from the end conditions o f the preceeding element. Finally, substituting the values of the nodal amplitudes A'~s and B ' s calculated from equations (25) and (26) into equations (15) and (16), the displacements and the rotati6ns o f any point along the b e a m at any time can be obtained.
Department of Civil Engineering, New Mexico State University, Las Cruces, New Mexico 88003, U.S.A.
O. L. ROUFAEIL
(Received 31 July 1985, and in revisedform 14 February 1986) REFERENCES
1. S. TIMOSHENKO 1921 Philosophical Magazine 41, 288-289. On the correction for shear of the differential equation for transverse vibration of prismatic bars. 2. R. W. TRAILL-NASH and A. R. COLLAR 1953 Quarterly Journal of Mechanics and Applied Mathematics 6, 186-222. The effect of shear flexibility and rotary inertia on bending vibration of beams. 3. T. C. HUANG 1961 Journal of Applied Mechanics 28, 579-584. The effect of rotary inertia and shear .deformation on the frequency of normal mode equations of beam with simple end conditions. 4. W. P. HOWSON and F. W. WILLIAMS 1973 Journal of Sound and Vibration 26, 503-515. Natural frequencies of frames with axially loaded Timoshenko members. 5. B. WENDORFF 1966 Theoretical Numerical Analysis. New York. Academic Press. 6. G. B. WARBURTON 1976 The Dynamical Behaviour of Structures. Oxford: Pergamon Press, second edition. 7. J. C. HOUBOLT 1950 Journal of Aeronautical Science 17, 540-550. A recurrence matrix solution for the dynamic response of elastic air craft. 8. K. J. BATHE 1982 Finite Element Procedure in Engineering Analysis. Englewood Cliffs, New Jersey: Prentice-Hall. 9. D. E. JOHNSON 1966 American Institute of Aeronautics and Astronautics Journal 4, 1450-1451. A proof of the stability of the Houbolt method.