Mathematics 174 – Applied Mathematics 2 Lab/Problem Set 1 – A “Shooting Party” January 27, 2000

Background

We have now seen that the Sturm-Liouville differential operator

L(ϕ) = \frac{d}{dx} (p(x)\frac{dϕ}{dx}) + q(x)ϕ,

restricted to the space of C^{2} functions ϕ satisfying the boundary conditions

\{ \begin{matrix}β_{1}ϕ(a) + β_{2}ϕ′(a) = 0 \\ β_{3}ϕ(b) + β_{4}ϕ′(b) = 0 \end{matrix} ,

is self-adjoint for the inner product

⟨f, g⟩ = ∫_{a}^{b} f(x)g(x) dx.

Recall that eigenfunctions are non-zero solutions of the equation L(ϕ) = −λσ(x)ϕ and the boundary conditions.

The self-adjointness implies that eigenfunctions for distinct eigenvalues are orthogonal with respect to the integral inner product with weight function σ(x). That is, if

\{ \begin{matrix}L(ϕ) = −λσ(x)ϕ(x) & and \\ L(ψ) = −μσ(x)ψ(x), & where \end{matrix}

λ ≠ μ, then we have

⟨ϕ, ψ⟩_{σ} = ∫_{a}^{b} ϕ(x)ψ(x)σ(x) dx = 0.

Also, assuming the Sturm-Liouville problem we consider is regular, there is a 1-dimensional space of eigenfunctions for each eigenvalue.

The other portions of the the Sturm-Liouville theorem (the form of the set of eigenvalues, the completeness of the eigenfunctions, and the pointwise convergence theorem for eigenfunction expansions) are more difficult to prove. We will look at a portion of the proof that the eigenvalues form a sequence λ_{n} → ∞ in the problems below, but we will not go into the details on the remaining portions.

Instead, today, we want to illustrate the theorem by looking at an example where the solutions of the Sturm-Liouville boundary value problem cannot be expressed in terms of the elementary functions from calculus. We will use numerical methods to graph approximate solutions and find the eigenvalues.

Lab Question

We will consider the Sturm-Liouville problem

\{ \begin{matrix}\frac{d^{2}ϕ}{dx^{2}} + λ(x^{2} + 1)ϕ = 0 \\ ϕ(0) = 0 \\ ϕ(1) = 0 \end{matrix}

where p(x) = 1, q(x) = 0, and σ(x) = x^{2} + 1. Since p(x) and σ(x) are continuous and strictly positive on [0, 1], this is a regular problem. Our goals are to approximate the first four eigenvalues λ_{1}, λ_{2}, λ_{3}, and λ_{4}, and to generate graphs of corresponding eigenfunctions ϕ_{1}, ϕ_{2}, ϕ_{3}, and ϕ_{4}. To do this we will use a rudimentary version of a numerical technique called the “shooting method”. The reason for the name should become clear as we proceed(!)

The basic idea of “shooting” is that we can solve the differential equation with initial conditions ϕ(0) = 0 and (for instance) ϕ′(0) = 1 for any given λ. Then, by adjusting λ, we can try to get the boundary condition ϕ(1) = 0 to hold as well. (Think of firing a cannon to hit a target at x = 1, where you can adjust the range using the value of λ!!)

1. For instance use the following Maple commands to load the DEtools package and define a function called shoot:

                          with(DEtools);
shoot := lambda -> DEplot(diff(phi(x),x$2) + lambda*(x^2+1)*phi(x),
      phi(x),x=0..1,[[phi(0)=0,D(phi)(0)=1]],stepsize=0.001);

Try entering the command shoot(5). The output is a plot of the solution of the differential equation with λ = 5. Then enter shoot(10). What do these graphs tell you about the location of the first eigenvalue: λ_{1}? Try zeroing in to find a good-looking approximation to λ_{1}.

2. Continue in the same way to generate approximate values of λ_{2}, λ_{3}, and λ_{4} together with graphs of eigenfunctions for each. Note the zeroes of the eigenfunctions in [0, 1].

3. Of course, the “eye-ball” method you used above is not terribly accurate. If we wanted better approximations to the eigenvalues, we would want to be able to examine the actual approximate value of the computed solution at x = 1 and see whether it is within some specified distance of 0 (e.g. |ϕ(1)| < 10^{−d} for some d). To do this in Maple, we can proceed as follows. First, we can define the solution of the initial value problem as a function (depending on λ, and also on x):

sol := lambda -> dsolve({diff(phi(x),x$2) + lambda*(x^2+1)*phi(x),
              phi(0)=0,D(phi)(0)=1},phi(x),numeric);

Then, for instance, phi := sol(5); defines the function ϕ which solves the differential equation and the initial conditions ϕ(0) = 0, ϕ′(0) = 1, with λ = 5. If you enter phi(1); the output

[x = 1, ϕ(x) = .2206662798196227, ϕ′(x) = −1.063068688690081]

gives the values of ϕ and its derivative at x = 1. From this we can see that λ = 5 is not an eigenvalue! Take your approximate eigenvalues from parts 1 and 2 and “polish” them a bit to get better approximations where the corresponding approximate eigenfunctions have |ϕ(1)| < 10^{−4}.

Problem Set Question

In this question, you will work through a big portion of a proof that a regular Sturm-Liouville boundary value problem has a sequence of real eigenvalues λ_{n} with λ_{n} → ∞ as n → ∞.

1) (Preparation – first order linear ODE). A first order linear ODE is a differential equation of the form

(1)  \frac{dy}{dx} = a(x)y + b(x)

a. The corresponding homogeneous equation is

(2)  \frac{dy}{dx} = a(x)y

show that the general solution of the homogeneous equation is

y_{h}(x) = y(0)e^{∫_{0}^{x} a(t) dt}.

b. Show that if y_{1}(x) and y_{2}(x) are any two solutions of (1), then y_{1}(x) − y_{2}(x) is a

solution of the homogeneous equation (2). c. Show that if we multiply both sides of

\frac{dy}{dx} − a(x)y = b(x)

by e^{− ∫ a(x) dx}, then the resulting equation can be integrated and solved for y to yield a particular solution of (1):

y_{p}(x) = e^{∫ a(x) dx} ∫ e^{− ∫ a(x) dx}b(x) dx

d. From parts a, b, c, deduce that the general solution of (1) can be written:

y(x) = e^{∫_{0}^{x} a(t) dt} (y(0) + ∫_{0}^{x} e^{− ∫_{0}^{s} a(t) dt}b(s) ds)

2) Consider the general (regular) Sturm-Liouville boundary value problem:

\{ \begin{matrix}L(ϕ) = \frac{d}{dx} (p(x)\frac{dϕ}{dx}) + q(x)ϕ = −λσ(x)ϕ \\ β_{1}ϕ(a) + β_{2}ϕ′(a) = 0 \\ β_{3}ϕ(b) + β_{4}ϕ′(b) = 0 \end{matrix}

Introduce notation X(x) = ϕ(x) and Y (x) = p(x)ϕ′(x). Show that the 2nd order differential equation can be rewritten as a system of two 1st order equations as follows:

\{ \begin{matrix}X′ = \frac{1}{p(x)}Y \\ Y′ = −(q(x) + λσ(x))X \end{matrix}

3) Make the following change to “polar-like” coordinates: X = R sin θ, Y = R cos θ. (Note: R and θ are corresponding functions of x.) Show that the system from part 2 becomes

\{ \begin{matrix}R′ sin θ + R cos θ · θ′ = \frac{1}{p}R cos θ \\ R′ cos θ − R sin θ · θ′ = −(q + λσ)R sin θ \end{matrix}

From this deduce the following two ODE’s:

\{ \begin{matrix}R′ = \frac{R}{2} (\frac{1}{p} − (q + λσ)) sin 2θ & (3) \\ θ′ = \frac{1}{p} cos^{2} θ + (q + λσ) sin^{2} θ & (4) \end{matrix}

4) Equation (4) does not depend on R. Hence general theory for differential equations implies that for each λ ∈ R (and any fixed c ∈ R), there is a solution θ_{λ}(x) of (4) with θ_{λ}(a) = c. Using this, show that the solution of (3) can be written as

R(x) = R(a)e^{∫_{a}^{x} F_{λ}(t) sin 2θ_{λ}(t) dt}

where F_{λ}(x) = \frac{1}{p(x)} − (q(x) + λσ(x)). Deduce that if R(a) ≠ 0, then X(x) = 0 only when θ_{λ}(x) = kπ for some integer k.

5) So, from part 4, it is enough to show that for each n ≥ 1, there is a λ_{n} such that θ_{λ_{n}}(b) = nπ. Those give the eigenvalues. (Note that this is essentially what were were doing in the lab portion with “shooting”.) The key step in showing this is the following Lemma:

Lemma. If λ < μ, and θ_{λ}(a) ≤ θ_{μ}(a), then θ_{λ}(x) < θ_{μ}(x) for all x ∈ (a, b].

To show this, let ω = θ_{μ}(x) − θ_{λ}(x). Using equation (4) above (which is satisfied by both θ_{μ} and θ_{λ}), find a 1st order linear ODE

ω′(x) = a(x)ω(x) + b(x)

satisfied by ω, then use part 1 above to determine the form of the solution with ω(a) = θ_{μ}(a) − θ_{λ}(a) > 0. Deduce the Lemma from the form of the solution of your ODE.

To complete the proof, there are still a few technical details to handle, but we will not consider them. The strategy is this:

i. First, by the Lemma, we know θ_{λ}(b) is increasing as a function of λ. ii. We would show that θ_{λ}(b) < π/2 for some λ, and that θ_{λ}(b) → ∞ as b → ∞. It

follows then that for each n, there is some λ_{n} such that θ_{λ_{n}}(b) = nπ.

Assignment

Individual writeups, due in class on Thursday, February 3.