This is an evaluated

(c) Copyright 1994,1995 by Evans M. Harrell, II. All rights reserved

assumptions in our derivation fo the wave equation are not valid. For example, the

constant c^2 = (spring constant)/(mass density) for the spring may be different at different

positions, leading to an equation of the form

2 2 2 2 2 \partial u/\partial t = (c[x]) \partial u/\partial x (ModWE1)

or, as a thorough examination of the derivation of the wave equation reveals, we could

more generally have:

2 2 \partial u/\partial t = p(x) (\partial /\partial x) s(x) \partial u/\partial x (ModWE2)for some potentially complicated positive functions p(x) and s(x). How well does the method of separation of variables do for problems like this? Rather well, actually, although we may have to encounter some new functions.

**Model Problem VI.3**.. Suppose that the wave speed depends on
position, so that

c^2 = 1/(1 + x), 0 < x < 1, DBC at 0 and 1.

Find the normal modes of vibration.

Solution.

When we attempt to solve the equation with the ansatz u[t,x] = T[t] X[x], we
find the

following eigenvalue problem.

- X''[x] = (1 + x) mu X[x]

There are actually some special functions, called Airy functions, which solve
the ODE

y''[x] = x y. Two independent solutions are called Ai(x) and Bi(x), or,
according to

*Mathematica* , AiryAi[x] and AiryBi[x].

*In[1]:=*

Plot[AiryAi[x], {x, -10,10}]

*In[2]:=*

Plot[AiryAi[x], {x, -5,3}]; Plot[AiryBi[x], {x, -5,3}]

*Out[2]=*

-Graphics-

The function Bi explodes exponentially to the right, while Ai decays
exponentially. They

both oscillate to the left (why?). By changing variables we can get these
functions to solve

our eigenvalue equation:

*In[3]:=*

D[AiryAi[-mu^(1/3) (x+1)],{x,2}]

*Out[3]=*

1/3 -(mu (1 + x) AiryAi[-(mu (1 + x))])

We need a linear combination Ai(- mu^(1/3) (1+x)) + C Bi ( - mu^(1/3) (1+x))
which is 0

at x=0 and 0 at x=1. I avoid the cube root at this stage by letting p = mu^1/3:

*In[4]:=*

FindRoot[{AiryAi[- p] + C AiryBi[- p] == 0, \ AiryAi[- 2 p] + C AiryBi[- 2 p] == 0}, \ {p,1.5},{C,2}]

*Out[4]=*

{p -> 1.87088, C -> 0.819688}

*In[5]:=*

Plot[AiryAi[-1.87088 (1+x)] +.819688 AiryBi[-1.87088 (1+x)], \ {x,0,1}]

*Out[5]=*

-Graphics-

The fact that this function has no nodes between 0 and 1, and therefore
resembles sin(9 x),

indicates that this is the spatial part of the fundamental (lowest-frequency)
mode. The

eigenvalue and normal mode are:

*In[6]:=*

mu0 = p /. %%; Mode0[t_,x_] = (A0 Cos[mu0 t] + B0 Cos[mu0 t]) *\ (AiryAi[-1.87088 (1+x)] +.819688 AiryBi[-1.87088 (1+x)])

*Out[6]=*

(AiryAi[-1.87088 (1 + x)] + 0.819688 AiryBi[-1.87088 (1 + x)]) (A0 Cos[1.87088 t] + B0 Cos[1.87088 t])

*In[7]:=*

FindRoot[{AiryAi[- p] + C AiryBi[- p] == 0, \ AiryAi[- 2 p] + C AiryBi[- 2 p] == 0}, \ {p,3},{C,-2}]