Journal of Sound and Vibration (1995) 180(4), 637–644
THE GENERATION OF WAVES IN INFINITE STRUCTURES BY MOVING HARMONIC LOADS M. S. A. H† Department of Engineering, University of Cambridge, Cambridge CB2 1PZ, England (Received 22 July 1993, and in final form 1 December 1993) The theory of convolution is extended to account for time-varying loads moving over infinite systems. Fourier transforms are used to simplify the convolution, reducing it to a multiplication of transforms of system impulse response and load. If a harmonic load is moving over the system it is found that the possible existence of travelling waves can be identified, for a given system, load frequency and velocity, without the need to perform the inverse Fourier transform, a task which is often difficult. The possible presence of travelling waves can be identified by a simple method involving straight line constructions on a plot of the system’s frequency spectrum. The phase velocities, group velocities and frequencies of waves ahead of and behind the load can be identified along with any critical speeds and velocities that may exist.
1. INTRODUCTION
Vibration waves in structures generated by moving or harmonic loads have been the subject of many investigations. Vibrations in beams, plates and shells were all studied by Lord Rayleigh [1, 2]. Timoshenko is credited with extending the theory for beams to allow for transverse shear [3], while the extension of simple plate theory for high frequency excitation was accomplished by Mindlin [4]. Vibration and wave propagation in layered media have been analyzed by Haskell [5] and more computer oriented techniques have been developed that can give practical solutions [6]. The response of a semi-infinite solid to a moving constant load was studied by Cole and Huth [7]. Further analysis of a generalized beam model on a Pasternak foundation [8] and generalized plate [9] has also been carried out for this loading. The analysis of such systems under constant moving loads has also been studied in detail and excellent reviews have been given elsewhere [10, 11]. Systems subjected to moving harmonic loads have also been studied in the literature. Mathews [12, 13] investigated the response of Euler beams on Winkler foundations. Chonan extended this to cover Timoshenko beam theory in the absence of damping [14]. Steady state motions of systems are almost universally sought by using Fourier transforms and inversion by using the residue theorem. This method is generalized here and a simple, graphical method is presented which identifies possible propagating modes of vibration and their presence ahead of or behind the moving load. While the magnitudes of the propagating modes are not identified by this method, it is possible to identify wave † Now at the Department of Civil Engineering and Building Science, University of Edinburgh, Edinburgh EH9 3JN, Scotland.
637 0022–460X/95/090637 + 08 $08.00/0
7 1995 Academic Press Limited
. . .
638
frequencies and speeds and critical load frequencies and speeds in a rapid way which also helps to interpret different regimes of system response. 2. MOVING CONVOLUTION THEORY
The response of a linear system to a time-varying input is given by the convolution integral [15] y(t) =
g
a
h(t − t) f (t) dt,
(1)
−a
where y(t) is the response at time t, f (t) is the input force at time t, and h(t) is the response at time t to a unit impulse at time t = 0. If the input is moving with respect to the system so that its position is given by the vector rf (t) at time t, and the response is measured at a point with position vector ry (see Figure 1), the impulse response function becomes dependent on three variables, and is given by h(ry , rf , t). The convolution integral in equation (1) then becomes y(ry , t) =
g
a
h(ry , rf (t), t − t) f (t) dt.
(2)
−a
This is a general expression for the response of a linear system to a moving input. It is frequently possible to make a number of simplifying assumptions, as follows. (i) If the system response is not a function of the output position but only of the separation of the input and output then h(ry , rf (t), t − t) = h(ry − rf (t), t − t).
(3)
(ii) If the input is moving at constant velocity then the origin may be chosen as the position of the input at time t = 0 and the problem is confined to one dimension: i.e., rf (t) = vt,
y(ry , t) = y(x, t),
h(ry − rf (t), t − t) = h(x − vt, t − t),
(4) (5)
where y(x, t) is the response at position x at time t, v is the vehicle speed, and h(x, t) is the response at position x and time t to a unit impulse at the origin at time t = 0.
Figure 1. Moving convolution integral definitions.
639
Under these conditions equation (2) becomes y(x, t) =
g
a
h(x − vt, t − t) f (t) dt.
(6)
−a
This equation is similar to that derived by Cebon [16], where it was used to find the response of a continuously supported beam to arbitrary, moving excitation. 3. FREQUENCY DOMAIN SOLUTION OF THE MOVING CONVOLUTION INTEGRAL
The simple convolution integral, equation (1), is often solved in the frequency domain by using Fourier transforms because the transform of the integral reduces to a simple multiplication [15]. The same technique can be used to simplify the moving convolution integral (equation (6)) but it is necessary to take Fourier transforms with respect to both time and space variables as follows: y˜ (j, v) =
=
0 1 gg 1 2p
2
a
y(x, t) e−ivt e−ijx dt dx
−a
0 1 ggg 1 2p
2
a
h(x − vt, t − t) f (t) e−ivt e−ijx dt dt dx
−a
= 2ph (j, v) f (v + vj).
(7)
Here v is the angular frequency of loading, corresponding to the time t, j is the wavenumber, corresponding to the distance x, and a tilde 0 indicates a transformed function. It is possible to predict, at least qualitatively, the behaviour of a system by examining the Fourier transform of its response. Equation (7) gives a simple method of doing this by examining the Fourier transforms of the system’s impulse response functions and the Fourier transform of the load separately and then combining them. If the load is varying harmonically, f (t) = f0 cos (v0 t), then its Fourier transform is simply f (v) = f0 (d(v − v0 ) + d(v + v0 )/2,
(8)
where d(v) is the Dirac delta function. This is most easily verified by applying the inverse transform. This can be substituted into equation (7) and the system response is given by the inverse Fourier transform y(x, t) = p
gg
a
h (j, v)(d(v + vj − v0 ) + d(v + vj + v0 )) eivt eijx jj dv.
(9)
−a
The Dirac delta functions make one of these integrals simple, leaving y(x, t) = p
g
a
(h (j, −vj + v0 ) ei(−vj + v0 )t + h (j, −vj − v0 ) ei(−vj − v0 )t ) eijx dj.
(10)
−a
If the Fourier transform of the impulse response function has poles at j = Jk and v = Vk(1) or Vk(2) , where Vk(1) = −vJk + v0 ,
and
Vk(2) = −vJk − v0 ,
(11)
. . .
640
then this integral may be evaluated by using the residue theorem. If the corresponding residues are Ak and Bk then the response is given by
0
y(x, t) = 2p 2i s Ak exp (i[Vk(1) t + Jk x]) + s Bk exp (i[Vk(2) t + Jk x])
0
k
k
1
1 1
= 2p 2i s Ak exp (i[v0 t + Jk (x − vt)] + s Bk exp (i[−v0 t + Jk (x − vt)] , (12) k
k
where the summations are taken over all poles above the real axis for x q vt and all poles below the real axis for x Q vt. If Jk is real then the corresponding exponential term represents a one-dimensional travelling wave with wavelength 2p/Jk and angular frequency cJk 2 v0 , where the negative sign corresponds to the term in the first summation and the positive sign to the term in the second summation. This method has been used to solve problems with moving, harmonic loads by many different authors. In general, it is not a simple matter to find the poles and residues of the Fourier transform of the impulse response function which satisfy equations (11). However, it is generally a much simpler matter to find combinations of real j and v for which h (j, v) is infinite and therefore for which a real Jk can exist. The combinations of j and v that are excited by the moving load are given by the arguments of the Dirac delta functions of equation (9) and travelling waves therefore exist if real solutions can simultaneously be found for j and v from the following equations: h (j, v) = a,
and either
v + vj + v0 = 0
or
v + vj − v0 = 0. (13a, b)
The travelling waves have frequency v and wavenumber j. From equation (12) these correspond to forward travelling waves if j and v have opposite signs. It is possible to interpret the simultaneous equations (13) graphically, and for systems such as layered strata a graphical solution is by far the simplest approach. If the load is not moving then the solutions to equations (13) are the usual, harmonically excited, travelling wave solutions. The set of points in (j, v) space for which these solutions exist is known as the frequency spectrum of the system [17] and is found from equation (13a) as this equation is concerned only with the system, not the loading. Examples of systems’ frequency spectra are given in the following sections and the literature mentioned in the introduction. The equations (13b) define straight lines in (j, v) space with gradients −v and intercepts on the v axis of 2v0 . These load lines are purely dependent on the velocity and frequency of the load. The intersections of the load lines and the spectrum give solutions of equations (13) and, therefore, are frequencies and wavenumbers of possible travelling waves. As j and v are required to have opposite signs for forward travelling waves it is convenient to plot the frequency spectrum and load lines on axes of j and v' = −v, as solutions in the positive quadrant of this space represent forward travelling waves and also the load line now has a gradient of +v.
4. THE RESPONSE OF AN EULER BEAM ON A WINKLER FOUNDATION
One of the simplest practical systems to which this theory may be applied is an infinite Euler beam resting on a Winkler foundation [16] as shown in Figure 2. It is necessary to investigate the impulse response of such a system.
641
Figure 2. An infinite beam on a Winkler foundation.
The equation of motion of a beam resting on a Winkler foundation, subjected to a point load moving over its surface is [10] EI 1 4w/1x 4 + m 1 2w/1t 2 + kw = d(x − d(t))P(t),
(14)
where w is the beam’s vertical displacement, d(x) is the Dirac delta function, d(t) is the position of the force at time t, EI is the beam’s bending stiffness, m is the beam’s mass per unit length and k is the foundation stiffness per unit length. The impulse response is obtained by setting d(t) = 0, and P(t) = d(t). Taking Fourier transforms with respect to both space and time variables gives the Fourier transform of the impulse response function h(x, t): h (j, v) = [(2p)2(EIj 4 − mv 2 + k)]−1 .
(15)
This is drawn on dimensionless axes in Figure 3 with a pair of typical load lines. There are four intersections between the load lines and frequency spectrum in the positive quadrant, marked A, B, C and D, indicating four possible travelling modes. The intersections marked A', B', C' and D' correspond to the same propagating modes as those in the positive quadrant but give amplitudes which are the complex conjugates. These solutions are required if equation (12) is to have a real solution. In order to decide which of these modes travel ahead of or behind the load it is necessary to consider the velocity of energy propagation associated with each mode, not the velocity of the propagation of peaks. Energy is propagated at the group velocity associated with each mode. The group velocity is given by the rate of change of frequency with wavenumber and is therefore given by the tangent to the frequency spectrum at the intersection of the frequency spectrum and load lines. The velocity of the load is given by the gradient of the load line and examination of the intersection quickly reveals which modes propagate energy more quickly than the load is moving, and therefore exist ahead of the load, and which modes propagate energy more slowly than the load is moving, and
Figure 3. The intersection of velocity lines (——) and frequency spectrum (
) for an Euler beam.
. . .
642
Figure 4. The presence of backward travelling waves.
therefore exist behind the load. With reference again to Figure 3, modes A and B propagate ahead of the load and modes C and D propagate behind it. All the modes propagate energy in the same direction as the moving load because dj/dv q 0. The phase velocity, i.e., the velocity of the peaks, of each mode is given by the gradient of a line joining the origin to the point of intersection. It is interesting that the phase velocity of mode B is slower than the velocity of the load so that it appears that the load is catching up with the waves that it generates. Similarly, mode D has a phase velocity that is faster than the velocity of the load so that it appears that these waves are catching up with the load. A stationary load is represented by a pair of horizontal lines. The lowest frequency for which waves will propagate (i.e., the cut-on frequency) is therefore vc = z(k/m). If the beam is excited at a frequency above the cut-on frequency then waves will be generated that travel away from the load in both directions. This is indicated by the intersections B and B' in the second and fourth quadrants of (j, v) space shown in Figure 4. In this illustration there are only two propagating modes present, A ahead of the load and B behind it. If the frequency of loading is kept constant and the speed increased it is possible to get four modes once again. Three of these modes will have positive group and phase velocities and one will have negative group and phase velocities. The critical speeds can be identified by drawing lines that are tangent to the spectrum and passing through the origin.
5. THE RESPONSE OF A TIMOSHENKO BEAM ON A WINKLER FOUNDATION
An Euler beam is not a good model to use at high frequencies and it is found that the Timoshenko beam models real behaviour more accurately [17]. The simultaneous partial differential equations describing the behaviour of a Timoshenko beam can be Fourier transformed to give the non-dimensional frequency equation V 4 − ((Q + R)X 2 + RQ + 1)V 2 + RQX 4 + RX 2 + RQ = 0,
(16)
where the non-dimensional groups are V = vzr/k ,
Q = kGA/zEIk,
4 X = jzEI/k,
R = (A/I)zEI/k,
Here r is the density of the beam material, k is the shear coefficient, G is the beam
643
Figure 5. The frequency spectrum for Timoshenko beam.
T 1 Timoshenko beam parameters Young’s modulus Shear modulus Second moment of area Cross-sectional area Density Timoshenko shear coefficient Winkler foundation stiffness
3000 MPa 1110 MPa 6·66 × 10−4 m4 0·2 m2 2300 kg/m3 0·833 235 MN/m2
E G I A r k k
material’s shear modulus and A is the beam’s cross-sectional area. All other parameters are as defined for the Euler beam above. Equation (16) is quadratic in V 2 and therefore has two positive solutions for each value of X giving rise to two complete curves as shown on Figure 5. The values of the parameters used to generate this frequency spectrum are given in Table 1 and are typical of values used to simulate asphalt road responses [18, 19]. The lower of the two curves in the frequency spectrum is due to the predominantly bending mode of wave propagation and the upper curve is due to the predominantly shear mode of wave propagation. The relevant cut-on frequencies are v1 = zk/rA
for bending,
and
v2 = zkAG/rI
for shear.
(17)
Construction of a load line on the frequency spectrum would indicate which modes are propagating and the tangent to the spectrum at the intersection would indicate if the mode propagates ahead of or behind the load, according to the rules outlined above.
6. CONCLUSIONS
The equations describing the frequency spectra of many systems may be found in the literature including those listed in the introduction. The method described is a rapid technique that can be used to investigate the existence of propagating waves causes by moving harmonic loads. It is also a valuable aid in detailed analysis as it indicates what regimes of behaviour are to be expected and the number of critical speeds/frequencies associated with a given system.
644
. . . REFERENCES
1. L R 1877 The Theory of Sound (two volumes). New York: Dover publications, second edition, 1945 re-issue. 2. L R 1887 Proceedings of the London Mathematical Society 17, 4–11. On waves propagated along the plane surface of an elastic solid. 3. S. P. T 1921 Philosophical Magazine 41, 744–746. On the correction for shear of the differential equation for transverse vibrations of prismatic bars. 4. R. D. M 1951 Journal of Applied Mechanics 18, 31–38. Influence of rotary inertia and shear on flexural motions of isotropic, elastic plates. 5. N. A. H 1953 Bulletin of the Seismological Society of America 43, 17–34. The dispersion of surface waves on multilayered media. 6. E. K and J. M. R¨ 1981 Bulletin of the Seismological Society of America, 17, 1743–1761. Stiffness matrices for layered soils. 7. J. C and J. H 1958 (December) Journal of Applied Mechanics, 433–436. Stresses produced in a half-plane by moving loads. 8. H. S and T. T 1980 Journal of Applied Mechanics 47, 879–883. Steady-state vibrations of a beam on a Pasternak foundation for moving loads. 9. M. A. E 1991 Journal of Physics D: Applied Physics 24, 541–546. Moving loads on an infinite strip of constant thickness. 10. L. F 1972 Vibration of Solids and Structures Under Moving Loads. Groningen, Netherlands: Noordhoff International. 11. A. D. K 1981 Solid Mechanics Archives 6(4), 401–449. Continuously supported beams and planes subjected to moving loads—a survey. 12. P. M. M 1958 Zeitschrift fu¨r Angewandte Mathematik und Mechanik 38(3/4), 105–115. Vibrations of a beam on elastic foundation I. 13. P. M. M 1959 Zeitschrift fu¨r Angewandte Mathematik und Mechanik 39(1/2), 13–19. Vibrations of a beam on elastic foundation II. 14. S. C 1978 Zeitschrift fu¨r Angewandte Mathematik und Mechanik 58(1), 9–15. Moving harmonic load on an elastically supported Timoshenko beam. 15. D. E. N 1985 Random Vibrations and Spectral Analysis London: Longman; second edition. 16. D. C 1988 Institution of Mechanical Engineers, Proceedings 202(C2), 109–117. Theoretical road damage due to dynamic tyre forces of heavy vehicles, part 2: simulated damage caused by a tandem-axle vehicle. 17. K. F. G 1975 Wave Motion in Elastic Solids. Oxford: Clarendon Press. 18. M. S. A. H and D. C 1992 Third International Symposium on Heavy Vehicle Weights and Dimensions, Cambridge, U.K., 65–72. Flexible pavement response models for assessing dynamic axle loads. London: Thomas Telford. 19. M. S. A. H 1990 Ph.D. Thesis, Cambridge University Engineering Department. The response of flexible pavements to dynamic tyre forces.