Unit 5: Boundary Value Problems for Partial Differential Equations
Chapter 24: Wave Equation
Sections24.5: longitudinal vibrations in an elastic rod
Copyright
Copyright * 2001 by Addison Wesley Longman, Inc.
All rights reserved. No part of this publication may be reproduced, stored in a retrieval system, or transmitted, in any form or by any means, electronic, mechanical, photocopying, recording, or otherwise, without the prior written permission of the publisher. Printed in the United States of America.
Intializations
> restart;
>
with(plots):
interface(showassumed=0);
Warning, the name changecoords has been redefined
>
Introduction
A rod of length
L
and density
(mass per unit volume) is deformed by pulling or pushing in opposite directions at its ends. The rod is said to stretch of compress
linearly
along its axis if the restorative forces induced in the rod obey Hooke's law. In effect, the rod behaves like a spring which is both extensible and compressible. Such a rod is called
elastic
if after a deformation, it returns to its original state.
Under compression, the distance between contiguous cross-Sections in the rod decreases, while under expansion (rerefaction), it increases. Compression or expansion experienced by neighboring plane cross-Sections can propagate through the rod as a wave which obeys the wave equation
Since the disturbance travels along the axis of the rod, the wave is said to be longitudinal . Thus, a steel rod, hammered on its end, will support longitudinal waves of compression and rarefaction traveling the length of the rod. Similarly, the rubber engine mounts supporting the engine in an automobile, when subjected to shock, experience the same passage of distortions progressing as plane waves through the medium.
Although the longitudinal vibrations in a rod are governed by the same wave equation as the vibrating string, the meaning of the displacement
is different, and subtly so. In the rod,
measures the displacement, at time
, of the plane in the rod that was at
initially. Since these displacements are relative to where the plane Sectionswas initially, there is an additivity which makes the complete motion much harder to fathom than the vibrations in the finite string where the displacement was directly observable.
In addition, the free endpoint condition whereby an end is free to move is formalized in terms of the strain , a new concept for the reader not familier with the theory of elasticity. Hence, our approach will be to state an appropriate boundary value problem, nurture intuition by illustrating the solution, and finally, develop the theory and derive the governing partial differential equation.
>
Example 24.6
Problem Statement
Formulate and solve an appropriate boundary value problem for the longitudinal elastic vibrations in a rod of length
whose left and right ends are both free. Assume
as the wave speed, and uniformly stretch the rod by one unit, releasing it without imparting any initial velocity.
>
Problem Formulation
If the rod is placed between
and
on an
-axis, and if a copy of the
-axis is etched on the rod at equilibrium, then the desired boundary value problem consists of the wave equation
the boundary conditions
and the initial conditions
=
The function
measures the displacement, at time
, of the plane Sectionsin the rod that was at location
initially. The endpoint conditions are derivatives because, as we will see later whenwe derive the governing equations, "free ends" mean no forces act on the ends. In fact, the elastic force on a face at location
x
is given by
.
To determine the initial distribution of displacements in the rod, we must first obtain the constant value of the strain , defined as the change in length, per unit length. Thus, we compute the constant strain as
=
The initial distribution of displacements is then found by accumulating the strain along the rod, obtaining
=
Finally, the initial velocity, namely,
, is zero.
This boundary value problem is new because we have not yet solved the wave equation with homogeneous Neumann conditions (as derivative conditions at the endpoints are called.)
>
Solution by Separation of Variables
Assume
, that is,
> U := X(x)*T(t);
>
as the separated form for the solution. As in Sections 24.1 and 24.2, the wave equation then yields the two ODEs
Apply this separation assumption to the endpoint conditions, obtaining
and
However, applying the separation assumption to the endpoint conditions now gives
=>
=>
=>
=>
In addition, applying the separation assumption to the initial velocity yields
=>
=>
a result we obtained each time we solved the wave equation with a similarly vanishing initial velocity.
We therefore face the Sturm-Liouville system
= 0
and the modified initial value problem
From Sections16.1, the Sturm-Liouville BVP for
has the eigenvalues
=
,
n
= 0, 1, 2, ...
and the corresponding eigenfunctions
,
n
= 0, 1, 2, ...
The corresponding solutions for
are then
,
n
= 0, 1, 2, ...
A typical eigensolution is therefore
and
is a linear combination of all such eigensolutions. Evaluating at
tells us the coefficients in the sum are the Fourier coefficients for
, so we have the solution
+
For
, that is, for
> f := x/Pi;
>
the Fourier coefficients are
and
which we obtain in Maple as follows. First, obtain
with the separate integral
> q0 := (2/Pi)*Int(f,x=0..Pi);
>
whose value is
> a0 := value(q0);
>
The general Fourier coefficient is given by
> qa := (2/Pi)*Int(f*cos(n*x),x=0..Pi);
>
Tell Maple that n is an integer via
> assume(n,integer);
>
and evaluate this defining integral, obtaining
> a := simplify(value(qa));
>
The partial sum
+ 2
that is,
> u20 := a0/2 + sum(a*cos(n*x)*cos(n*t),n=1..20):
>
is graphed as a surface in Figure 24.18, below.
> plot3d(u20, x=0..Pi, t=0..4*Pi, axes=frame, color=black, labels=[x,t,u], labelfont=[TIMES,ITALIC,12], scaling=constrained, tickmarks=[3,6,2], orientation=[-45,50]);
>
The wave-like properties of the solution are apparent, but the detailed interpretation of the meaning of these waves requires more investigation. Remember,
gives the displacement, at time
, of the plane Sectionsthat was at
initially. What we want to see is a representation of the physical motion of the elastic rod, and unlike the transverse waves on the finite string,
does not directly show these physical motions.
An animation of the solution
is a dynamic interpretation of the solution surface, but is still not a representation of the physical motions in the rod. It shows the displacements
which each face that was at location
is now experiencing at time
. It tells us that whatever is happening is happening linearly, and it takes place from the endpoints in, and back, with a wave-front like a shock wave. Several frames from this animation appear in Figure 24.19, but the full animation is generated in Maple by executing the following command.
>
p1 := animate(u20, x=0...Pi, t=0..4*Pi, frames=80, labels=[x,u], color=black, labelfont=[TIMES,ITALIC,12], xtickmarks=3, ytickmarks=[0,1], scaling=constrained):
p1;
>
What does this animation mean? What is it telling us about the longitudinal vibrations in the rod? For openers, the animation shows the displacements
which each face that was at location
x
is now experiencing at time
t
. Unfortunately, it is still hard from this animation to comprehend the actual motions taking place in the rod.
However, it does tell us that whatever is happening is happening linearly, and it takes place from the endpoints in, and back, with a wave-front like a shock wave.
A more revealing animation simulates the actual rod by imagining a set of eleven equispaced scratch marks on the rod, and following the motion of these scratch marks in real time. Figure 24.20 shows the bar at
but the full animation is generated by the following Maple commands.
>
uu := unapply(u20, x, t):
F1 := (X,T) -> [[X+uu(X,T),0],[X+uu(X,T),1]]:
F2 := T -> {seq(F1(k*Pi/10,T), k = 0..10)}:
F3 := T -> plot(F2(T), color=black):
p2 := display([seq(F3(k*Pi/20), k = 0..79)], insequence=true, scaling=constrained, ytickmarks=2, labels=[x,``]):
p2;
>
The uniformly stretched rod relaxes to its natural length, with the collapse occurring from the "outside, in." The center of the rod is the only part of the rod that does not move. Then, the collapse continues past the point of reaching natural length, and the rod contracts by an amount equal to the original stretch. This continuing contraction now takes place from "inside, out" as the center of the rod first experiences the compression, with the outside of the rod feeling the compression last. The rod then expands from this state of full compression to natural length, from "outside, in." Having achieved natural length, the rod continues expanding to its initial stretched state, this time from the "inside, out."
In addition to this simulation of the rod itself, it is also instructive to examine an animation of
, which is the localized (per-face) value of the strain,
. This animation shows more clearly how the compression and rarefaction waves travel with a "front" through the rod.
To compute the derivative
we have to differentiate the infinite sum representing the solution
. The derivative of a Fourier series converges even more slowly than the original series. Where we obtained acceptable results with a twenty-term approximation for
itself, it will take 100 terms for the derived series to converge. This translates into a modest increase in computation time.
Figure 24.21 shows, at times
and
, the strain in the rod. However, executing the following Maple commands will generate the complete animation.
>
u100 := a0/2 + sum(a*cos(n*x)*cos(n*t),n=1..100):
ux100 := diff(u100,x):
p3 := animate(ux100, x=0..Pi, t=0..4*Pi, frames=80, color=black, labels=[x,strain], xtickmarks=3, ytickmarks=[-.3,0,.3]):
p3;
>
The animation shows, as a function of time, the changing strain in the rod. The strain moves longitudinally through the rod, as a shock wave. It starts out uniform across the rod, and drops to zero with a shock-front. When the strain is completely gone, the rod has reached its natural length. As the rod then enters a compressive mode, the strain becomes negative, and the region of compression expands from the center out to the ends. At this juncture the strain graph is a constant negative amount.
The compression then starts to vanish from the ends inwards until the rod again reaches its natural length. This corresponds to the vanishing of the strain. As the rod enters an expansion phase, the strain appears above the x -axis, growing, with a shock-front, outwards from the center. The expansion continues until the rod returns to its original stretched state.
Finally, we place side-by-side the three animations of displacement, physical motion, and strain. Run continuously, and slowly enough, the resulting animation allows the three facets of the solution to be compared and coordinated.
> display(array([p||(1..3)]));
>
Theory and Derivation
A rod which obeys a linear law of elasticity is said to obey
Hooke's law
, wherein forces are proportional to displacement. This linear relationship between force and displacement, expressed by the formula
, is a staple of elementary calculus and general physics, especially in any discussion of the behavior of a spring. In fact, the linear spring was an essential element in Chapters 5 and 6. Thus,
, the force required to stretch a spring, is proportional to
, the amount of the stretch. The constant of proportionality,
, is called the spring constant.
In the context of elasticity, Hooke's law is sometimes stated as
stress is proportional to the strain
where stress is the force per unit cross-sectional area and
strain
=
is defined as the elongation per unit length. In the context of elasticity, Hooke's law is then
stress = E strain
where E is Young's modulus, the constant of proportionality between stress and strain . This gives Young's modulus the units of
E ~
~
=
In fact, starting with
, and dividing through by
A
, the cross-sectional area, we have
The displacement
is the elongation
from the definition of strain, namely,
. Hence, the spring law becomes
=
=
=
where
, Young's modulus, is
. Hence, the two statements of Hooke's law are equivalent.
Since longitudinal waves travel in the direction of the rod's axis, we want a coordinate system along that axis. However, we must assume a linear elasticity which keeps plane Sections both plane and parallel. Thus, forces within the rod obey the linear Hooke's law and preserve the orientation of all plane Sections which are perpendicular to the axis of the rod.
Put the rod (whose equilibrium length is
) on a ruler, aligning the left end of the rod with the left end of the ruler. At equilibrium, paint a copy of the ruler's coordinates on the rod. For clarity, take this linear coodinate system to be an
x
-coordinate system.
Identify a face (i.e., a plane section) within the rod by the coordinate it had at equilibrium. If the equilibrium coordinates are painted on the rod, then these coordinates travel with the rod should the rod be stretched or compressed. In addition, leave the ruler on the table so the original meaning of an x -coordinate is preserved.
Consider
R
, a segment of the rod bounded at equilibrium by faces at
x
and
. This segment has length
and is shown in the following figure, FIgure 24.22.
>
p4 := plot([10+2*cos(t),3+3*sin(t),t=0..2*Pi],color=black):
p5 := plot([2*cos(t),3+3*sin(t),t=Pi/2..3*Pi/2],color=black):
p6 := plot([3+2*cos(t),3+3*sin(t),t=Pi/2..3*Pi/2],color=black):
p7 := plot([5+2*cos(t),3+3*sin(t),t=Pi/2..3*Pi/2],color=black):
p8 := plot([5+2*cos(t),3+3*sin(t),t=-Pi/2..Pi],color=black, linestyle=2):
p9 := plot(6,x=0..10,color=black):
p10 := plot(0,x=0..12,color=black):
p11 := textplot({[2,-.5,`x`], [5,-.5,`x+dx`], [12.5,0,`x`], [10,-.5,`L`]}, font=[TIMES,ITALIC,12]):
p12 := textplot([2.4,3,`R`], font=[TIMES,ITALIC,12]):
display([p||(4..12)],scaling=constrained,view=[-2..13,-1..6],axes=none);
>
Let
denote the displacement experienced, at time
t
, by the face, which at equilibrium, was located at
x
. It helps to remember that this face still bears the coordinate
x
. The displacement undergone by this face is the quantity
, and is measured along the axis of the rod. Unlike the vibrating string where the displacements are very obvious, the lateral displacements of faces in the rod are much more difficult to visualize.
The partial differential equation governing longitudinal vibrations in the rod will be a consequence of
F
=
m
a
, Newton's second law of motion. It's not hard to guess that
a
, the acceleration of the face with the coordinate
x
painted on it, is given by
. Determining
F
, the net force acting on this face, is more subtle.
At time t , the left and right faces of the segment R are now located respectively at
>
X[left] := x+u(x,t);
X[right] := (x+dx)+u(x+dx,t);
>
The length of the deformed segment is then
> Delta := X[right]-X[left];
>
The change in length of this segment is therefore the relative displacement between the left and right faces of the segment, namely
> Delta - dx;
>
the relative displacement between the left and right faces of the segment. The strain associated with the face which, at equilibrium, was at location x , is given by the limiting ratio
Notice that strain, an extensive property, is ascribed to a point by this limit. This definition parallels the definition of pressure which is also an extensive property, force per unit area, ascribed to a point.
Thus, the strain on a single face with coordinate x painted on it is
> q3 := Limit((u(x+dx,t)-u(x,t))/dx,dx=0);
>
so that the value of the strain is
, that is,
> strain = convert(value(q3),diff);
>
This is a fundamental result in the theory of elasticity. Getting this right at the beginning is well worth the time invested in getting it. Thus, at time
and location
, the strain caused by the displacement
, is positive if the applied force tends to increase the size of the region R, and negative otherwise.
The force on a plane Sectionscan be obtained from the strain, using Hooke's law, stress = E strain . With A as the cross-sectional area of the rod, we have
=
stress
=
E
strain
so that
force
=
E
A
strain
=
E
A
Thus, the elastic force on the face whose equilibrium coordinate was x , is given by
To represent this force in Maple, define
> F := (x,t) -> E*A*D[1](u)(x,t);
>
so that
> 'F(x,t)' = convert(F(x,t),diff);
>
The net force on segment R is given by the sum of the forces acting on the left and right faces of R , namely,
which, is Maple, is
> `net force` := F(x+dx,t)-F(x,t);
>
The total mass of segment R is given by
since the volume of the segment is A , the cross-sectional area, times the length, dx .
Taking the acceleration of a face as
, we can write the
i
-component of Newton's second law,
F
=
m
a
, as
that is, as
> q4 := `net force` = rho*A*dx*diff(u(x,t),t,t);
>
Divide by
to get
that is,
> q5 := simplify(q4/(A*rho*dx));
>
which, in the limit as
, becomes
that is,
> q6 := convert(map(limit,q5,dx=0),diff);
>
Upon defining
, the partial differential equation becomes
, that is
> algsubs(E/rho=c^2,q6);
>
With a slight rearrangement and the introduction of subscript notation for derivatives, this equation can be written as
the very same wave equation we derived for the vibrating string. The difference here is that
represents the displacement of the face that was at
x
at time
t
, and this displacement is lateral, along the axis of the rod.
>