Module 3: Series Solution of Differential Equations

Chapter 1: Introduction to Series Solutions

Learning Objectives
By the end of this chapter, learners will be able to:

● Understand the concept of power series and its role in solving differential equations.
● Define radius of convergence and apply the ratio test.
● Apply power series methods to solve second-order differential equations.
● Construct and solve recurrence relations arising from substitution.
● Interpret general solutions in series form.

1.1 Introduction to Series Solutions

Elementary functions (such as exponentials, trigonometric, logarithmic, or polynomial functions)
cannot be used to solve many differential equations. In these situations, we look for a series
solution, which means we attempt to represent the answer as an infinite sum of the independent
variable's powers. Assuming that a differential equation's solution may be expressed as a power
series and then substituting it into the equation to find the unknown coefficients is a potent
strategy.

Second-order linear differential equations of the following type, where P(x) and Q(x) are analytic
functions close to a point x = x₀ , can be solved using this method especially well:
d²y/dx² + P(x) dy/dx + Q(x) y = 0

1.2 What is a Power Series?

An infinite sum of the following type is called a power series:
From n = 0 to ∞,
y(x) = Σ [aₙ * (x − x₀ )ⁿ]
where x₀ is the series' centre, which is frequently assumed to be 0,
aₙ are constants (coefficients),
and x is the variable.

Though it may be limitless in length, this is comparable to a polynomial. An illustration of a

power series with x = 0 at its centre:
y(x) is equal to a₀ + a₁ x + a₂ x² + a₃ x³ +...

1.3 Radius of Convergence

For all values of x, a power series does not necessarily converge. The interval in which the
series converges is defined by the radius of convergence, or R.
If the series is:
Σ aₙ (x − x₀ )ⁿ

then the series converges absolutely when:

|x − x₀ | < R

and diverges when:

|x − x₀ | > R

The radius R can often be found using the ratio test:

R = lim (as n → ∞) |aₙ / aₙ₊ ₁ |

Within this radius, we can use the power series as a valid solution to the differential equation.

1.4 Strategy for Power Series Solutions

Suppose a second-order ODE has been given of the form:
y'' + P(x) y' + Q(x) y = 0
and a solution is to be found near x = 0.
The general procedure for solving this:

1. Assume a power series solution:

y(x) = Σ (from n = 0 to ∞) [aₙ xⁿ]

2. Differentiate term-by-term to find y' and y'':

y'(x) = Σ (from n = 1 to ∞) [n * aₙ * xⁿ⁻ ¹]
y''(x) = Σ (from n = 2 to ∞) [n * (n − 1) * aₙ * xⁿ⁻ ²]

3. Substitute these expressions into the original ODE.

4. Rearrange the resulting expression as a single power series.

5. Equating coefficients of like powers of x to zero yields a recurrence relation among the
coefficients aₙ.

6. Solve the recurrence to obtain the coefficients (a₀ , a₁ , ...), and thus the power series
solution.

1.5 Example: Solving y'' − x·y = 0

Example Equation:
y'' − x y = 0

We aim to find a solution around x = 0 using a power series.

Step 1: Assume a Power Series Solution

Let:
y(x) = Σ (from n = 0 to ∞) [aₙ xⁿ]
y'(x) = Σ (from n = 1 to ∞) [n * aₙ xⁿ⁻ ¹]
 y''(x) = Σ (from n = 2 to ∞) [n * (n − 1) * aₙ xⁿ⁻ ²]

Step 2: Plug into the Differential Equation

The equation is:
y'' − x y = 0

Substitute the expressions:

Σ (from n = 2 to ∞) [n(n − 1) * aₙ xⁿ⁻ ²] − x * Σ (from n = 0 to ∞) [aₙ xⁿ] = 0

Simplify the second term:

x * Σ aₙ xⁿ = Σ (from n = 0 to ∞) [aₙ xⁿ⁺ ¹]

So the equation becomes:

Σ (from n = 2 to ∞) [n(n − 1) * aₙ xⁿ⁻ ²] − Σ (from n = 0 to ∞) [aₙ xⁿ⁺ ¹] = 0

Step 3: Reindex Both Series to a Common Power of x

In the first term, let k = n − 2

So when n = 2, k = 0
Then n = k + 2
Thus:
Σ (from n = 2 to ∞) [n(n − 1) * aₙ xⁿ⁻ ²]
= Σ (from k = 0 to ∞) [(k + 2)(k + 1) * aₖ ₊₂ xᵏ]

Second term:
Σ (from n = 0 to ∞) [aₙ xⁿ⁺ ¹]
= Σ (from k = 1 to ∞) [aₖ ₋₁ xᵏ]
(Shift index: let k = n + 1 ⇒ n = k − 1)

Now combine both series:

Σ (from k = 0 to ∞) [(k + 2)(k + 1) * aₖ ₊₂ xᵏ] − Σ (from k = 1 to ∞) [aₖ ₋₁ xᵏ] = 0

Handle the k = 0 term separately:

At k = 0:
[(0 + 2)(0 + 1) * a₂ x⁰] = 2 * a₂

Then for k ≥ 1:
[(k + 2)(k + 1) * aₖ₊₂ − aₖ₋₁] xᵏ

2 * a₂ + Σ (from k = 1 to ∞) {[(k + 2)(k + 1) * aₖ₊₂ − aₖ₋₁] xᵏ} = 0

Step 4: Equating Coefficients of xᵏ to Zero

Since this series equals zero for all x, each coefficient must be zero:
For k = 0:
2 * a₂ = 0 ⟹a₂ = 0

For k ≥ 1:
 (k + 2)(k + 1) * aₙ₊ ₂ − aₙ₋ ₁ = 0
⟹aₖ ₊₂ = aₖ ₋₁ / [(k + 2)(k + 1)]

This is the recurrence relation.

Step 5: Use Recurrence to Build the Series

Assume a₀ and a₁ are arbitrary constants corresponding to initial conditions

Already found:
a₂ = 0

Now use recurrence:

For k = 1:
a₃ = a₀ / [3 * 2] = a₀ / 6
For k = 2:
a₄ = a₁ / [4 * 3] = a₁ / 12
For k = 3:
a₅ = a₂ / [5 * 4] = 0
For k = 4:
a₆ = a₃ / [6 * 5] = (a₀ / 6) / 30 = a₀ / 180
For k = 5:
a₇ = a₄ / [7 * 6] = (a₁ / 12) / 42 = a₁ / 504

Thus:
y(x) = a₀ [1 + (1/6) x³ + (1/180) x⁶ + ...] + a₁ [x + (1/12) x⁴ + (1/504) x⁷ + ...]

This is the general solution expressed as a power series

1.6 Check Your Progress

Q1. A power series centered at has the general form:
A) ∑ₙ₌ ₀ ^∞ cₙ · (x + a)ⁿ
B) ∑ₙ₌ ₀ ^∞ cₙ · (x − a)ⁿ
C) ∑ₙ₌ ₀ ^∞ (xⁿ − aⁿ)
D) ∑ₙ₌ ₁ ^∞ (x − a)ⁿ⁻ ¹

Q2. If a function y(x) is analytic at x = a, it means:

A) y(x) has a derivative only at x = a
B) y(x) is defined only on a finite interval
C) y(x) can be represented by a power series around x = a
D) y(x) cannot be represented by a power series at all

Q3. The interval in which a power series converges is known as:

A) Interval of accuracy
B) Interval of smoothness
C) Interval of convergence
D) Radius of expansion

Q4. The radius of convergence R of a power series

∑ₙ₌ ₀ ^∞ cₙ · (x − a)ⁿ
can be found using:
A) limₙ→∞ |cₙ|
B) limₙ→∞ 1⁄|cₙ|
C) limₙ→∞ |cₙ⁄cₙ₊ ₁ |
D) ∑ |cₙ|

Q5. Which of the following is not a valid solution method for second-order linear differential
equations?
A) Variation of parameters
B) Power series method
C) Frobenius method
D) Euler’s Forward method

Q6. The power series solution to a differential equation assumes the solution:
A) Cannot be expressed in any closed form
B) Is a finite polynomial
C) Can be represented as an infinite series of powers of x
D) Is a constant function

Q7. In the method of power series, to find coefficients, we typically:

A) Use matrices
B) Use Laplace transforms
C) Differentiate term-by-term and substitute into the DE
D) Integrate both sides of the DE

Q8. Which condition must be satisfied to use the power series method at a point x = a?
A) The point must be a singular point
B) The function must be discontinuous at x = a
C) The differential equation must be linear and analytic at x = a
D) The point must be a turning point

Q9. The power series

∑ₙ₌ ₀ ^∞ (xⁿ⁄n!)
is the Taylor expansion of:
A) ln(1 + x)
B) sin x
C) eˣ
D) tan⁻ ¹ x

Q10. If a power series solution to a differential equation is

y(x) = ∑ₙ₌ ₀ ^∞ aₙ · xⁿ,
then the value of a₀ is usually determined by:
A) Initial condition y′(0)
B) Boundary condition
C) Initial condition y(0)
D) Maximum value of the function

Answer Key
1. B
2. C
3. C
4. C
5. D
6. C
7. C
8. C
9. C
10. C

1.7 Solved Illustrations

Illustration 1

Problem:
Solve the differential equation
d²y/dx² − y = 0
using the power series method about x = 0.

Solution:

Assume a solution of the form:

y(x) = ∑ₙ₌ ₀ ^∞ aₙ · xⁿ
=> y(x) = a₀ + a₁ x + a₂ x² + a₃ x³ + ...
Differentiate term-by-term:

dy/dx = ∑ₙ₌ ₁ ^∞ n · aₙ · xⁿ⁻ ¹

d²y/dx² = ∑ₙ₌ ₂ ^∞ n · (n − 1) · aₙ · xⁿ⁻ ²

Now substitute into the differential equation:

d²y/dx² − y = 0
=> ∑ₙ₌ ₂ ^∞ n(n − 1)aₙ · xⁿ⁻ ² − ∑ₙ₌ ₀ ^∞ aₙ · xⁿ = 0

Change index of the first sum so that both sums are powers of xⁿ:

Let m = n − 2
Then n = m + 2
So:

∑ₙ₌₂^∞ n(n − 1)aₙ · xⁿ⁻²

= ∑ₘ₌₀^∞ (m + 2)(m + 1)aₘ₊₂ · xᵐ
Now both sums become:

∑ₘ₌₀^∞ [(m + 2)(m + 1)aₘ₊₂ − aₘ] · xᵐ = 0

Since this is 0 for all x, equating coefficients:

(m + 2)(m + 1)aₙ₊ ₂ − aₙ = 0
=> aₙ₊ ₂ = aₙ⁄[(m + 2)(m + 1)]

Now compute terms using initial coefficients a₀ and a₁ :

a₂ = a₀ ⁄2
a₃ = a₁ ⁄6
a₄ = a₂ ⁄(4·3) = (a₀ ⁄2)⁄12 = a₀ ⁄24
a₅ = a₃ ⁄(5·4) = (a₁ ⁄6)⁄20 = a₁ ⁄120

Hence the power series becomes:

y(x) = a₀ + a₁ x + (a₀ ⁄2)x² + (a₁ ⁄6)x³ + (a₀ ⁄24)x⁴ + (a₁ ⁄120)x⁵ + ...

Group even and odd terms:

y(x) = a₀ (1 + x²⁄2! + x⁴⁄4! + ...) + a₁ (x + x³⁄3! + x⁵⁄5! + ...)

Which is:

y(x) = a₀ · cosh(x) + a₁ · sinh(x)

Final Answer: y(x) = a₀ · cosh(x) + a₁ · sinh(x)

Illustration 2
Problem:
Find the power series solution for the differential equation
(1 − x²) · d²y/dx² − x · dy/dx + y = 0
about x = 0 up to 4 non-zero terms.

Solution:

This is Legendre’s Equation with n = 0. Let’s proceed with:

Assume y(x) = ∑ₙ₌ ₀ ^∞ aₙ · xⁿ

=> dy/dx = ∑ₙ₌ ₁ ^∞ n · aₙ · xⁿ⁻ ¹
=> d²y/dx² = ∑ₙ₌ ₂ ^∞ n(n − 1) · aₙ · xⁿ⁻ ²

Substitute into the given equation:

(1 − x²) · d²y/dx² − x · dy/dx + y = 0

Compute each term:

● (1 − x²) · d²y/dx² = ∑ₙ₌ ₂ ^∞ n(n − 1) · aₙ · xⁿ⁻ ² − ∑ₙ₌ ₂ ^∞ n(n − 1) · aₙ · xⁿ

● −x · dy/dx = −∑ₙ₌ ₁ ^∞ n · aₙ · xⁿ
● y = ∑ₙ₌ ₀ ^∞ aₙ · xⁿ

Combine and rearrange terms as a single power series:

We obtain:
∑ₙ₌ ₀ ^∞ [ (n + 2)(n + 1)aₙ₊ ₂ − n(n − 1)aₙ + n · aₙ − aₙ ] · xⁿ = 0

Equating coefficients of like powers of xⁿ:

(n + 2)(n + 1) · aₙ₊ ₂ − n(n − 1) · aₙ − n · aₙ + aₙ = 0

=> (n + 2)(n + 1) · aₙ₊ ₂ = [n(n − 1) + n − 1] · aₙ
=> (n + 2)(n + 1) · aₙ₊ ₂ = (n² − 1) · aₙ
=> aₙ₊ ₂ = [(n² − 1)⁄((n + 2)(n + 1))] · aₙ

Let’s compute the terms assuming a₀ = 1, a₁ = 0 (for even function):

a₂ = [(0² − 1)⁄(2·1)] · a₀ = (−1⁄2) · 1 = −1⁄2

a₄ = [(2² − 1)⁄(4·3)] · a₂ = (3⁄12) · (−1⁄2) = −1⁄8

So the series becomes:

y(x) = a₀ + a₂ x² + a₄ x⁴ + ...
=> y(x) = 1 − (1⁄2)x² − (1⁄8)x⁴ + ...

Final Answer:
y(x) ≈ 1 − (1⁄2)x² − (1⁄8)x⁴
Illustration 3

Problem:
Solve the differential equation
dy/dx = y,
with the initial condition y(0) = 1
using the power series method.

Solution:

We assume the solution is of the form:

y(x) = ∑ (n = 0 to ∞) aₙ xⁿ
= a₀ + a₁ x + a₂ x² + a₃ x³ + ⋯

Differentiating term-by-term:
dy/dx = ∑ (n = 1 to ∞) n aₙ xⁿ⁻ ¹
= a₁ + 2a₂ x + 3a₃ x² + ⋯

But dy/dx = y(x)

So, equating both power series:

a₁ + 2a₂ x + 3a₃ x² + ⋯ = a₀ + a₁ x + a₂ x² + ⋯

Comparing coefficients of like powers of x:

x⁰: a₁ = a₀
x¹: 2a₂ = a₁
x²: 3a₃ = a₂
x³: 4a₄ = a₃
⋯

Using y(0) = a₀ = 1:

Then,
a₁ = a₀ = 1
a₂ = a₁ / 2 = 1/2
a₃ = a₂ / 3 = 1/6
a₄ = a₃ / 4 = 1/24
⋯
aₙ = 1 / n!

So the power series becomes:

y(x) = ∑ (n = 0 to ∞) xⁿ / n!
= eˣ

Answer: y(x) = eˣ

Illustration 4
Problem:
Solve the differential equation
y'' - y = 0
with initial conditions y(0) = 1 and y'(0) = 0
using power series.

Solution:

Assume y(x) = ∑ (n = 0 to ∞) aₙ xⁿ
Then,
y''(x) = ∑ (n = 2 to ∞) n(n - 1) aₙ xⁿ⁻ ²

So,
y''(x) - y(x) = ∑ (n = 2 to ∞) n(n - 1) aₙ xⁿ⁻ ² - ∑ (n = 0 to ∞) aₙ xⁿ

We re-index the first series:

Let m = n - 2 ⇒ n = m + 2

Then,
∑ (m = 0 to ∞) (m + 2)(m + 1) a₍ₘ₊₂₎ xᵐ - ∑ (n = 0 to ∞) aₙ xⁿ = 0

Combine both series:

∑ (n = 0 to ∞) [(n + 2)(n + 1) a₍ ₙ₊ ₂ ₎ - aₙ] xⁿ = 0

For the series to vanish identically, each coefficient must be zero:

(n + 2)(n + 1) a₍ ₙ₊ ₂ ₎ - aₙ = 0
⇒ a₍ ₙ₊ ₂ ₎ = aₙ / [(n + 2)(n + 1)]

Now use initial conditions:

a₀ = y(0) = 1
a₁ = y'(0) = 0

Compute other coefficients:

a₂ = a₀ / (2×1) = 1/2
a₃ = a₁ / (3×2) = 0
a₄ = a₂ / (4×3) = (1/2) / 12 = 1/24
a₅ =0
a₆ = a₄ / (6×5) = (1/24) / 30 = 1/720
⋯

So,
y(x) = a₀ + a₁ x + a₂ x² + a₃ x³ + a₄ x⁴ + a₅ x⁵ + a₆ x⁶ + ⋯
= 1 + 0 + (1/2)x² + 0 + (1/24)x⁴ + 0 + (1/720)x⁶ + ⋯
= ∑ (n = 0 to ∞) x²ⁿ / (2n)!

This is the Maclaurin series for cosh(x).

Answer: y(x) = cosh(x)

Illustration 5

Problem:
Solve y'' + y = 0 using a power series centered at x = 0, with y(0) = 0 and y'(0) = 1.

Solution:

Assume y(x) = ∑ (n = 0 to ∞) aₙ xⁿ
Then,
y''(x) = ∑ (n = 2 to ∞) n(n - 1) aₙ xⁿ⁻ ²

So:
y''(x) + y(x) = 0
⇒ ∑ (n = 2 to ∞) n(n - 1) aₙ xⁿ⁻ ² + ∑ (n = 0 to ∞) aₙ xⁿ = 0

Re-index the first term: Let m = n - 2 ⇒n = m + 2

⇒∑ (m = 0 to ∞) (m + 2)(m + 1) a₍ₘ₊₂₎ xᵐ + ∑ (n = 0 to ∞) aₙ xⁿ = 0
⇒∑ (n = 0 to ∞) [(n + 2)(n + 1) a₍ₙ ₊₂₎ + aₙ ] xⁿ = 0

⇒ (n + 2)(n + 1) a₍ ₙ₊ ₂ ₎ + aₙ = 0
⇒ a₍ ₙ₊ ₂ ₎ = -aₙ / [(n + 2)(n + 1)]

Initial conditions:
y(0) = a₀ = 0
y'(0) = a₁ = 1

Then:
a₂ = -a₀ / (2×1) = 0
a₃ = -a₁ / (3×2) = -1/6
a₄ = -a₂ / (4×3) = 0
a₅ = -a₃ / (5×4) = 1/120
⋯

So,
y(x) = a₀ + a₁ x + a₂ x² + a₃ x³ + a₄ x⁴ + ⋯
= 0 + x + 0 - (1/6)x³ + 0 + (1/120)x⁵ + ⋯
= ∑ (n = 0 to ∞) [(-1)ⁿ x²ⁿ⁺ ¹ / (2n + 1)!]

This is the Maclaurin series for sin(x).

Answer: y(x) = sin(x)

1.8 Unsolved Problem Set

1. Solve the following differential equation using the power series method centered at x =
0:

(i) y'' + y = 0

2. Find a power series solution centered at x = 0 for the equation:

(ii) y'' − x y = 0

3. Determine the first four non-zero terms of the power series solution about x = 0 for:

(iii) x² y'' + y = 0

4. Solve the following second-order differential equation using power series around x = 0:

(iv) (1 − x²) y'' − 2x y' + 2y = 0

5. Find the general solution in the form of a power series about x = 0 for the equation:

(v) y'' + x y' + y = 0

6. Given the differential equation y'' − 2x y' + 4y = 0, find the recurrence relation and the
first four terms of the power series solution about x = 0.

7. Use the method of power series to solve:

(vii) y'' + (1 + x²) y = 0

8. Solve using a power series centered at x = 0:

(viii) y'' − y' − x y = 0

9. Find the power series solution of:

(ix) y'' − x² y = 0

10. Consider the equation y'' + x² y = 0. Find the recurrence relation and compute the first
five non-zero terms of the solution.
Chapter 2: Theory of Power Series Methods

Learning Objectives
After studying this chapter, learners will be able to:

● Understand the theoretical basis of the power series method for solving second-order
differential equations.
● Differentiate between ordinary and singular points.
● Apply the existence and uniqueness theorem for power series solutions.
● Compute the radius of convergence using standard tests.
● Use analytic continuation to extend series solutions beyond their initial interval of
convergence.
● Recognize the limitations of power series methods.

2.1 Overview
One of the most basic methods for solving differential equations is the Power Series Method,
particularly in situations when it is difficult to find solutions in terms of simple functions. For
linear second-order ordinary differential equations (ODEs), this approach is quite effective.
Understanding the structure and convergence behaviour of power series, as well as the
circumstances in which a power series solution exists and is unique, is crucial for applying the
power series approach systematically. By establishing pertinent theorems, the radius of
convergence, and behaviour close to ordinary and singular points, this chapter develops the
theoretical foundation for the approach.

2.2 Power Series Review

A power series is an infinite series of the form:

y(x) = ∑ₙ₌ ₀ ^∞ aₙ (x - x₀ )ⁿ

Here:

● aₙ are the coefficients of the series

● x₀ is the center of the series
● x is the variable

This series converges within a radius R around x₀, known as the radius of convergence.

The radius R is determined by the Cauchy–Hadamard formula:

1 / R = lim supₙ→∞ |aₙ|^(1/n)

Or in practical terms, one can often use the Ratio Test:

 R = limₙ→∞ |aₙ / aₙ₊ ₁ |

If this limit exists and is non-zero, then the series converges absolutely for all |x - x₀| < R.

2.3 Application to Differential Equations

Suppose we are given a second-order linear differential equation of the form:

y'' + P(x) y' + Q(x) y = 0

We aim to find a power series solution of the form:

y(x) = ∑ₙ₌ ₀ ^∞ aₙ (x - x₀ )ⁿ

The method proceeds as follows:

1. Compute y'(x) and y''(x) from the assumed series form.

2. Substitute y, y', and y'' into the given differential equation.
3. Align powers of (x - x₀) and combine like terms.
4. Set the coefficient of each power of (x - x₀) equal to zero. This leads to a recurrence
relation.
5. Use the recurrence relation to find successive coefficients aₙ.

2.4 Ordinary and Singular Points

A point x = x₀ is called an ordinary point of the differential equation if the functions P(x) and
Q(x) are analytic at x = x₀, i.e., can be expanded as a power series around x₀.

If either P(x) or Q(x) is not analytic at x = x₀, then the point is called a singular point.

Definition:

● An ordinary point is where a power series solution is guaranteed to exist.

● A singular point may still admit a solution, but different methods such as Frobenius
method are used.

2.5 Existence and Uniqueness Theorem

Theorem:
If P(x) and Q(x) are analytic at a point x = x₀, then the differential equation:
 y'' + P(x) y' + Q(x) y = 0

has a unique power series solution of the form:

y(x) = ∑ₙ₌ ₀ ^∞ aₙ (x - x₀ )ⁿ

that converges at least within some radius R > 0 around x₀.

This theorem guarantees that a convergent power series solution around an ordinary point is
always found.

2.6 Radius of Convergence of Series Solution

The radius of convergence of the power series solution to a differential equation is at most the
distance from the center x₀ to the nearest singular point of the equation.

Let x₁ be the closest singular point to x₀, then:

R = |x₁ - x₀ |

This is based on the fact that the power series solution fails to converge at or beyond a singular
point.

2.7 Strategy for Solving Second-Order ODEs via Power Series

Consider a second-order linear ordinary differential equation of the form:

y'' + P(x) y' + Q(x) y = 0

Assume that the functions P(x) and Q(x) are analytic at the point x = x₀ .

Then, x = x₀ is called an ordinary point.

We attempt a solution in the form of a power series centered at x = x₀ :

y(x) = Σ (n = 0 to ∞) aₙ (x - x₀ )ⁿ

Differentiating term-by-term:

y'(x) = Σ (n = 1 to ∞) n·aₙ (x - x₀ )ⁿ⁻ ¹

y''(x) = Σ (n = 2 to ∞) n·(n - 1)·aₙ (x - x₀ )ⁿ⁻ ²

These expressions for y, y', and y'' are substituted into the original differential equation. The
resulting infinite series is simplified, and coefficients of like powers of (x - x₀ ) are equated to
zero.

This gives a recurrence relation for the coefficients aₙ, which can be used to construct the full
series solution.

2.8 Illustrative Example

Solve the differential equation:

y'' - x·y = 0

Assume a power series solution centered at x = 0:

y(x) = Σ (n = 0 to ∞) aₙ xⁿ

Then, compute the derivatives:

y''(x) = Σ (n = 2 to ∞) n·(n - 1)·aₙ xⁿ⁻ ²

Also,

x·y(x) = x · Σ (n = 0 to ∞) aₙ xⁿ = Σ (n = 0 to ∞) aₙ xⁿ⁺ ¹

Now substitute into the equation y'' - x·y = 0:

Σ (n = 2 to ∞) n·(n - 1)·aₙ xⁿ⁻ ² - Σ (n = 0 to ∞) aₙ xⁿ⁺ ¹ = 0

To combine these, we reindex the summations to write both series as powers of xᵏ.

Reindexing the Series

First term:

Let k = n - 2,

so n = k + 2

Then:

Σ (k = 0 to ∞) (k + 2)·(k + 1)·aₖ ₊₂ xᵏ

Second term:

Let k = n + 1,

so n = k - 1
Then:

Σ (k = 1 to ∞) -aₖ ₋₁ xᵏ

Now, rewrite the equation as:

Σ (k = 0 to ∞) (k + 2)(k + 1)aₖ ₊₂ xᵏ - Σ (k = 1 to ∞) aₖ ₋₁ xᵏ = 0

We want to write this as a single sum:

Let us split the first term at k = 0:

At k = 0: The first term becomes (2)(1)a₂ x⁰ = 2a₂

So we now write:

[2a₂] + Σ (k = 1 to ∞) [(k + 2)(k + 1)aₖ₊₂ - aₖ₋₁] xᵏ = 0

Equating Coefficients of Powers of x

For the equation to be true for all x, each coefficient of xᵏ must be zero.

So:

● For k = 0:

2a₂ = 0

a₂ = 0

● For k ≥ 1:

(k + 2)(k + 1)aₙ₊ ₂ - aₙ₋ ₁ = 0

This gives the recurrence relation:

aₙ₊ ₂ = aₙ₋ ₁ / [(k + 2)(k + 1)]

This recurrence formula allows generation of all further terms aₙ once the values of a₀ and a₁
are known.

Understanding the Structure of the Solution

now construct the terms of the series step-by-step using the recurrence relation.

Assume:

a₀ = A (arbitrary constant)
a₁ = B (arbitrary constant)
It is known that:

a₂ = 0

Now compute the rest:

● For k = 1:
a₃ = a₀ / [(3)(2)] = A / 6

● For k = 2:
a₄ = a₁ / [(4)(3)] = B / 12

● For k = 3:
a₅ = a₂ / [(5)(4)] = 0 / 20 = 0

● For k = 4:
a₆ = a₃ / [(6)(5)] = (A / 6) / 30 = A / 180

And so on.

Thus, the series looks like:

y(x) = A[1 + (1/6)x³ + (1/180)x⁶ + ...] + B[x + (1/12)x⁴ + ...]

This is the general solution expressed as a power series in x.

2.9 Radius of Convergence

When we solve a differential equation using a power series expansion, we need to ensure that
the resulting series actually converges. The interval within which it converges is governed by the
radius of convergence (R).

Let a solution be expressed as:

y(x) = Σ (from n=0 to ∞) aₙ (x - x₀ )ⁿ

The radius of convergence is calculated using the Ratio Test:

R = lim (as n → ∞) |aₙ / aₙ₊ ₁ |

The series converges for values of x such that:

|x - x₀ | < R

And diverges for:

|x - x₀ | > R
At |x - x₀| = R, the convergence behavior needs to be tested separately.

Note: The radius of convergence is often determined by the nearest singular point (i.e., where
the differential equation ceases to be analytic).

2.10 Analytic Continuation

Sometimes it is necessary to know how the answer behaves outside of a certain range, yet the
solution obtained from a power series is only valid inside that range.

Analytic continuation is a technique for expanding a power series' domain by creating new
series centred at various locations that coincide with the initial interval of convergence. It is
useful in achieving this.

For instance, if the original series is centered at x₀ and converges in the range |x - x₀| < R, we
can construct another series centered at x₁ such that:

|x₁ - x₀ | < R

and then extend from there.

This technique is particularly helpful in solving second-order differential equations near singular
points, as traditional series expansions typically break down.

2.11 Limitations of Power Series Method

While the power series method is a powerful tool, it has certain limitations:

1. Requires coefficients to be analytic functions:

The method assumes that P(x), Q(x), and R(x) (from the differential equation
P(x)y'' + Q(x)y' + R(x)y = 0) are analytic at the point of expansion.

2. Breaks down at singular points:

If x₀ is a singular point of the equation, convergence may not to seen in the series. In
such cases, the Frobenius method is employed.

3. Convergence may be limited:

If a series converges, it is possible that such convergence may only be seen within a
small radius, limiting practical usability.

4. No general solution format:

Many cases lack a closed-form expression for the sum. This requires truncation at a
finite number of terms for approximation.
2.12 Key Terms
● Power Series
● Radius of Convergence
● Ordinary Point
● Singular Point
● Recurrence Relation
● Analytic Continuation

2.13 Summary
● Power series offer a systematic method for solving second-order ODEs.
● Solutions exist and are unique near ordinary points.
● The radius of convergence is limited by the nearest singularity.
● Analytic continuation helps expand the validity of a solution.

2.14 Check Your Progress

1. The general form of a linear second-order differential equation to which the power series
method can be applied is:

A) y′ + P(x)y = Q(x)
B) y′′ + y′ + y = 0
C) P(x)y″ + Q(x)y′ + R(x)y = 0
D) y″ = f(x)

2. A function f(x) can be represented by a power series about a point x₀ if and only if:

A) f(x) is continuous at x₀
B) f(x) is differentiable at x₀
C) f(x) is integrable on an interval containing x₀
D) f(x) is analytic at x₀

3. Which of the following statements is true regarding the radius of convergence of a power
series?

A) It depends only on the degree of the differential equation

B) It is infinite for all power series
C) It equals the distance to the nearest singularity from the expansion point
D) It equals the coefficient of the highest power of x
4. The power series
Σ (from n=0 to ∞) aₙ (x − x₀)ⁿ
is centered at:

A) x = 0
B) x = 1
C) x = x₀
D) x = aₙ

5. The power series method is best suited when the coefficients of the differential equation are:

A) Constant
B) Polynomial
C) Non-continuous
D) Analytic at the point of expansion

6. If y(x) = Σ (from n=0 to ∞) aₙ xⁿ is a solution of a differential equation, then the term-by-term

differentiation is valid when:

A) The function is even

B) The series has only a finite number of non-zero terms
C) The series converges uniformly in an interval
D) The interval is bounded

7. The ratio test for convergence of a power series involves calculating:

A) lim (as n→∞) |aₙ₊ ₁ − aₙ|

B) lim (as n→∞) |aₙ₊ ₁ / aₙ|
C) lim (as n→∞) aₙ
D) lim (as n→∞) aₙ / n

8. For the series solution to be valid, the point x = x₀ must be:

A) A singular point
B) A regular point
C) A periodic point
D) A critical point

9. The power series expansion method is generally not applicable when:

A) The point of expansion is regular
B) The equation is second order
C) The equation contains singularities at the point of expansion
D) The function has a Maclaurin expansion

10. The interval of convergence is defined as:

A) The maximum value of the power series

B) The region where the series converges absolutely
C) The sum of the power series
D) The domain of the function

Answer Key:
1. C
2. D
3. C
4. C
5. D
6. C
7. B
8. B
9. C
10. B

2.15 Solved Illustrations

Illustration 1: Determining the Radius of Convergence

Problem:

Determine the radius of convergence of the power series

S(x) = ∑ (n = 1 to ∞) [ (3ⁿ * xⁿ) / n ].

Solution:

use the root test or Cauchy-Hadamard formula:

Let aₙ = (3ⁿ / n), then the radius of convergence is given by:

R = 1 / lim sup |aₙ|^(1/n)

= 1 / lim sup [(3ⁿ / n)^(1/n)]
= 1 / lim sup [3 * (1/n)^(1/n)]
As n → ∞, (1/n)^(1/n) → 1
So, R = 1 / 3

Answer: Radius of Convergence is 1/3

Illustration 2: Power Series Solution of a Differential Equation

Problem:

Find the first four terms of the power series solution for the differential equation:
y'' − y = 0

Assume a solution of the form y(x) = ∑ (n = 0 to ∞) aₙ * xⁿ

Solution:

Start by assuming:

y(x) = ∑ aₙ * xⁿ

y''(x) = ∑ (n = 2 to ∞) n(n − 1) * aₙ * xⁿ⁻ ²

Substitute into the equation:

∑ (n = 2 to ∞) n(n − 1) * aₙ * xⁿ⁻ ² − ∑ (n = 0 to ∞) aₙ * xⁿ = 0

Shift indices in first summation to match powers:

Let m = n − 2

n=m+2

So, first sum becomes:

∑ (m = 0 to ∞) (m + 2)(m + 1) * aₘ₊₂ * xᵐ

Now both sums are:

∑ (m = 0 to ∞) [ (m + 2)(m + 1) * aₘ₊₂ − aₘ ] * xᵐ = 0

Equating coefficients of like powers to zero:

(m + 2)(m + 1) * aₙ₊ ₂ − aₙ = 0
aₙ₊ ₂ = aₙ / (m + 2)(m + 1)

Initial terms: Let a₀ = A, a₁ = B (arbitrary constants)

Then: a₂ = a₀ / 2 = A / 2
a₃ = a₁ / 6 = B / 6
a₄ = a₂ / 12 = (A / 2) / 12 = A / 24
a₅ = a₃ / 20 = (B / 6) / 20 = B / 120

So the series becomes:

y(x) = A * [1 + (1/2)x² + (1/24)x⁴ + ...] + B * [x + (1/6)x³ + (1/120)x⁵ + ...]

Answer: y(x) = A * ∑ (n = 0 to ∞) [x²ⁿ / (2n)!] + B * ∑ (n = 0 to ∞) [x²ⁿ⁺ ¹ / (2n + 1)!]

Which is the series expansion of: y(x) = A * cosh(x) + B * sinh(x)

Illustration 3: Power Series Solution of y'' + y = 0

Problem:
Find the power series solution of the differential equation
y'' + y = 0
with initial conditions y(0) = 1, y'(0) = 0.

Solution:

Assume

y(x) = ∑ (n = 0 to ∞) aₙ * xⁿ

y''(x) = ∑ (n = 2 to ∞) n(n − 1) * aₙ * xⁿ⁻ ²

Substitute into the differential equation:

∑ (n = 2 to ∞) n(n − 1) * aₙ * xⁿ⁻ ² + ∑ (n = 0 to ∞) aₙ * xⁿ = 0

Shift index of the first summation:

Let m = n − 2

n=m+2

Then:
∑ (m = 0 to ∞) (m + 2)(m + 1) * aₘ₊₂ * xᵐ + ∑ (m = 0 to ∞) aₘ * xᵐ = 0

Combine the sums:

∑ (m = 0 to ∞) [ (m + 2)(m + 1) * aₘ₊₂ + aₘ ] * xᵐ = 0

Equating coefficients of xᵐ:

(m + 2)(m + 1) * aₙ₊ ₂ + aₙ = 0

aₙ₊ ₂ = −aₙ / (m + 2)(m + 1)

Start with a₀ = 1 (given), a₁ = 0 (given)

Then: a₂ = −a₀ / (2 * 1) = −1 / 2
a₃ = −a₁ / (3 * 2) = 0
a₄ = −a₂ / (4 * 3) = (1 / 2) / 12 = 1 / 24
a₅ = −a₃ / (5 * 4) = 0
a₆ = −a₄ / (6 * 5) = −(1 / 24) / 30 = −1 / 720

So,

y(x) = a₀ + a₁ x + a₂ x² + a₃ x³ + a₄ x⁴ + a₅ x⁵ + a₆ x⁶ + ...
= 1 − (1/2)x² + (1/24)x⁴ − (1/720)x⁶ + ...

This is the Maclaurin series for cos(x).

Final Answer:
y(x) = cos(x)

Illustration 4: Finding Radius and Interval of Convergence

Problem:
Determine the radius and interval of convergence of the series
S(x) = ∑ (n = 1 to ∞) [ (−1)ⁿ * xⁿ / n² ]

Solution:

Let aₙ = (−1)ⁿ / n²
Then:
S(x) = ∑ aₙ * xⁿ

We apply the ratio test:

L = lim (n → ∞) |aₙ₊ ₁ * xⁿ⁺ ¹ / (aₙ * xⁿ)|

= lim (n → ∞) |(−1)ⁿ⁺ ¹ / (n + 1)² * xⁿ⁺ ¹ ÷ ((−1)ⁿ / n² * xⁿ)|
= lim (n → ∞) |x * n² / (n + 1)²|

Now:
n² / (n + 1)² → 1 as n → ∞
L = |x| * 1 = |x|

The series converges when L < 1 ⇒ |x| < 1

Thus, radius of convergence R = 1

Now test at the endpoints:

1. At x = 1:
S(1) = ∑ (−1)ⁿ / n²
 Converges by alternating series test

2. At x = −1:
S(−1) = ∑ (−1)ⁿ * (−1)ⁿ / n² = ∑ 1 / n²

Converges

Final Answer:
Radius of convergence R = 1
Interval of convergence: [−1, 1]

2.16 Unsolved Problem Set

Section A: Power Series Solutions Around Ordinary Points

1. Solve the following differential equation using a power series about x = 0:
y'' − x * y = 0

2. Find the power series solution for:

y'' + 2x * y' + 3y = 0, with initial conditions y(0) = 1, y'(0) = 0

3. Obtain the first four non-zero terms of the power series solution to:
y'' − 3x * y' + 2y = 0, given y(0) = 2, y'(0) = 1

4. Determine the power series expansion of the solution to:

y'' + y = 0, about x = π/2

5. Solve the differential equation using the method of power series:

(1 − x²) * y'' − 2x * y' + 2y = 0

Section B: Radius and Interval of Convergence

6. Find the radius and interval of convergence of the power series:
∑ (n = 1 to ∞) [ xⁿ / (n * 3ⁿ) ]

7. For the series ∑ (n = 0 to ∞) [ (−1)ⁿ * xⁿ / n! ], determine its radius of convergence and

whether it converges uniformly on ℝ.

8. Determine the radius of convergence of the series:

∑ (n = 1 to ∞) [ x²ⁿ / n² ]
 9. Find the radius of convergence of the series solution to:
y'' − y = 0, centered at x = 1

10. Given the series ∑ (n = 1 to ∞) [ (x − 2)ⁿ / n³ ], determine its interval of convergence.

Section C: Method of Frobenius

11. Use the Frobenius method to find the solution near x = 0 for:
x² * y'' + x * y' − y = 0

12. Apply Frobenius method to solve:

x² * y'' + 2x * y' + x² * y = 0

13. Solve the differential equation using Frobenius method:

x² * y'' + (1 − x) * y' + y = 0

14. Find the indicial equation and roots for:

x² * y'' + x * y' + (x − 1) * y = 0

15. Solve using Frobenius method:

x² * y'' + x * y' + (x² − ν²) * y = 0 (Bessel’s Equation)

Section D: Application and Conceptual Reasoning

16. Show that the power series method fails at a singular point and explain why Frobenius is
necessary.

17. Explain the difference between regular and irregular singular points using examples.

18. For the equation x² * y'' + 5x * y' + 3y = 0, classify the point x = 0 and determine an
appropriate method of solution.

19. Derive the recurrence relation for the general solution of:
y'' − x * y' + y = 0, and identify the pattern.

Prove that the power series solution to y'' = y leads to the exponential function, and determine
the radius of convergence..

References
● Boyce, W. E., & DiPrima, R. C. (2017). Elementary Differential Equations and Boundary
Value Problems (11th ed.). Wiley.
● Simmons, G. F. (1972). Differential Equations with Applications and Historical Notes.
McGraw-Hill.
● Coddington, E. A. (1955). An Introduction to Ordinary Differential Equations. Prentice
Hall.
Chapter 3: Frobenius Method, Legendre Equation, and Bessel
Functions

Learning Objectives
After completing this chapter, learners will be able to:

● Understand the Frobenius method and its application to differential equations with
regular singular points.
● Derive and solve the indicial equation.
● Distinguish between different cases based on the roots of the indicial equation.
● Apply the Frobenius method to obtain power series solutions.
● Solve examples involving Frobenius expansions.

3.1 Introduction to Frobenius Method

The Frobenius Method is an extension of the power series method used to find solutions to
second-order linear differential equations near a regular singular point. It is generally used when
the power series method is incapable of generating a solution or provides a solution limited to a
particular range.

While the power series method assumes a solution of the form:

y(x) = Σ (n = 0 to ∞) aₙ·xⁿ

the Frobenius method further generalizes this to:

y(x) = Σ (n = 0 to ∞) aₙ·xⁿ⁺ ʳ
where r is a number to be determined (possibly non-integer).

This approach is essential when dealing with differential equations having a singularity at the
point of expansion.

3.2 Standard Form for Frobenius Method

Consider the linear second-order differential equation:

x²·y'' + x·p(x)·y' + q(x)·y = 0

If x = 0 is a regular singular point, then p(x) and q(x) are analytic at x = 0, or have poles of order
1 and 2 respectively.

The Frobenius method seeks a solution of the form:

 y(x) = xʳ·Σ (n = 0 to ∞) aₙ·xⁿ = Σ (n = 0 to ∞) aₙ·xⁿ⁺ ʳ

This series is substituted into the differential equation and to collect terms to determine the
coefficients aₙ and the exponent r.

3.3 Indicial Equation

The lowest power of x in the substituted expression provides the indicial equation. It determines
the possible values of r.

This equation is obtained by collecting terms of the lowest power of x after substitution:

r(r - 1) + p₀ ·r + q₀ = 0

This is a quadratic equation in r and is called the indicial equation.

Here,

p₀ = lim (x → 0) x·p(x) and q₀ = lim (x → 0) x²·q(x).

The roots r₁ and r₂ of the indicial equation determine the form of the solution. The nature of
these roots (real and distinct, equal, or differing by an integer) affects the strategy used for
obtaining a second linearly independent solution.

3.4 Cases Based on Indicial Roots

Let r₁ and r₂ be the roots of the indicial equation with r₁ ≥ r₂ .

Case 1: Distinct Roots, Not Differing by an Integer (r₁ - r₂ ∉ ℕ)

● Two linearly independent solutions exist of the form:

y₁ (x) = xʳ¹·Σ (n = 0 to ∞) aₙ·xⁿ y₂ (x) = xʳ²·Σ (n = 0 to ∞) bₙ·xⁿ

Case 2: Distinct Roots Differing by an Integer (r₁ - r₂ ∈ ℕ)

● The smaller root may yield a solution that contains log(x) terms.
● The first solution:
y₁ (x) = xʳ¹·Σ (n = 0 to ∞) aₙ·xⁿ
● The second solution is of the form:
y₂ (x) = C·y₁ (x)·ln(x) + xʳ²·Σ (n = 0 to ∞) bₙ·xⁿ

Case 3: Equal Roots (r₁ = r₂ )

● One solution is as above. The second solution usually involves log(x).

3.5 Steps of Frobenius Method
1. Rewrite the differential equation in standard form.
2. Assume a solution of the form

y(x) = Σ aₙ·xⁿ⁺ ʳ.

3. Differentiate and substitute y, y', and y'' into the equation.

4. Collect powers of x and set the coefficient of the lowest power to zero to solve for the
indicial equation.
5. Solve the indicial equation for r
6. Recursively solve for aₙ using the recurrence relation

3.6 Example: Solving via Frobenius Method

Consider the equation:

x²·y'' + x·y' - y = 0

This is a Cauchy-Euler equation with a regular singularity at x = 0.

Assume:

y = Σ aₙ·xⁿ⁺ ʳ

Then:

y' = Σ (n = 0 to ∞) (n + r)·aₙ·xⁿ⁺ ʳ⁻ ¹

y'' = Σ (n = 0 to ∞) (n + r)(n + r - 1)·aₙ·xⁿ⁺ ʳ⁻ ²

Substitute into the equation and simplify to get the indicial equation:

r² - 1 = 0

r = ±1

Then proceed to find aₙ using recurrence.

Key Terms
● Frobenius Method
● Regular Singular Point
● Indicial Equation
● Recurrence Relation
● Logarithmic Solution
Summary
● Frobenius method generalizes the power series to handle singular points.
● Assumes solutions of the form y(x) = xʳ·Σ aₙ·xⁿ
● Indicial equation determines allowable values of r.
● Nature of r₁ and r₂ defines solution structure.
● Used when ordinary power series fails or is insufficient.

Check Your Progress

Q1. Which of the following is the correct form of a Frobenius series solution about a regular
singular point ?

A. y(x) = ∑ (n = 0 to ∞) aₙ xⁿ
B. y(x) = ∑ (n = 0 to ∞) aₙ xⁿ⁻ ¹
C. y(x) = ∑ (n = 0 to ∞) aₙ xⁿ⁺ ʳ
D. y(x) = ∑ (n = 1 to ∞) aₙ xⁿ

Q2. The Frobenius method is typically used when the point under consideration is:

A. An ordinary point
B. A singular point of order 2
C. A regular singular point
D. An irregular singular point

Q3. What is the indicial equation in the Frobenius method?

A. A differential equation in aₙ
B. A quadratic equation in x
C. A polynomial equation in r
D. A recurrence relation in n

Q4. Suppose the indicial equation has two roots and such that is not an integer. Then:

A. Only one solution exists

B. The second solution is not guaranteed
C. Two linearly independent Frobenius solutions exist
D. Solution must be found numerically
Q5. If the indicial roots and differ by a positive integer, the second solution may involve:

A. Cosine and sine terms

B. Bessel functions
C. A logarithmic term
D. Legendre polynomials

Q6. In the Frobenius method, the coefficients are determined using:

A. Matrix algebra
B. A recurrence relation
C. Euler’s method
D. Separation of variables

Q7. Which of the following is not an assumption made in the Frobenius method?

A. The point is a regular singular point

B. The solution is in the form of a power series
C. The solution has only integer powers of x
D. The lowest power of x is not necessarily an integer

Q8. The indicial equation is obtained by:

A. Solving the ODE directly

B. Substituting the series into the ODE and equating the lowest power of x
C. Using Euler’s theorem
D. Using variation of parameters

Q9. For the equation

x² y″ + x y′ − y = 0,
identify the indicial equation when solving by Frobenius method at .

A. r(r − 1) = 0
B. r² = 0
C. r² − 1 = 0
D. r² − r = 0

Q10. Let the indicial roots be , . Then the general solution (if both series converge) is of the
form:

A. y(x) = A x² + B x² log(x)
B. y(x) = A x² + B x¹
C. y(x) = A x² + B x² sin(x)
D. y(x) = A x² + B x⁴

Answer Key

1. C
2. C
3. C
4. C
5. C
6. B
7. C
8. B
9. D
10. B

Solved Illustrations

Illustration 1:

Problem:
Solve using Frobenius method:
x²·y'' + x·y' - y = 0

Step 1: Identify standard form and regular singular point

Given: x²·y'' + x·y' - y = 0

This is already in standard form:

x²·y'' + x·p(x)·y' + q(x)·y = 0

Here,

p(x) = 1 ⇒ x·p(x) = x,

so p₀ = lim (x → 0) x = 0
q(x) = -1

x²·q(x) = -x²,

so q₀ = lim (x → 0) -x² = 0

So, the indicial equation becomes:

r(r - 1) + p₀ ·r + q₀ = 0
Step 2: Find the indicial equation
r(r - 1) = 0

r = 0 or r = 1

solve for both roots.

Step 3: General series form

Assume a solution of the form:

y(x) = Σ (n = 0 to ∞) aₙ·xⁿ⁺ ʳ

Then:

y' = Σ (n = 0 to ∞) (n + r)·aₙ·xⁿ⁺ ʳ⁻ ¹
y'' = Σ (n = 0 to ∞) (n + r)(n + r - 1)·aₙ·xⁿ⁺ ʳ⁻ ²

Step 4: Substitute into equation

substitute into the original ODE:

x²·y'' + x·y' - y = 0 becomes:

Σ (n = 0 to ∞) (n + r)(n + r - 1)·aₙ·xⁿ⁺ ʳ

+ Σ (n = 0 to ∞) (n + r)·aₙ·xⁿ⁺ ʳ

- Σ (n = 0 to ∞) aₙ·xⁿ⁺ ʳ = 0

Combining:

Σ (n = 0 to ∞) [(n + r)(n + r - 1) + (n + r) - 1]·aₙ·xⁿ⁺ ʳ = 0

So,

[(n + r)(n + r - 1) + (n + r) - 1] = (n + r)² - 1

Hence,

Σ (n = 0 to ∞) [(n + r)² - 1]·aₙ·xⁿ⁺ ʳ = 0

Step 5: Recurrence relation

We now get:

[(n + r)² - 1]·aₙ = 0

To satisfy the above for all n, we require:

For a₀ ≠ 0, the coefficient must vanish at n = 0, giving the indicial equation:

(r)² - 1 = 0

r = ±1

Let us solve for r = 1.

Then recurrence becomes:

[(n + 1)² - 1]·aₙ = 0 ⇒ (n² + 2n + 1 - 1)·aₙ = n(n + 2)·aₙ = 0

So for n ≠ 0,
aₙ = 0 for all n ≥ 1

Hence,

y₁ (x) = a₀ ·x¹ = a₀ ·x

Now try r = -1

Then,

[(n - 1)² - 1]·aₙ = (n² - 2n)·aₙ = n(n - 2)·aₙ = 0

So nonzero coefficients only for:

● n = 0 ⇒ a₀
● n = 2 ⇒ a₂ (since for n = 1, a₁ = 0)

Set a₀ = arbitrary, and find a₂ .

We substitute into the recurrence relation:

For n = 2:

(2 - 1)² - 1 = 1 - 1 = 0

a₂ is free (arbitrary)

So second solution:

y₂ (x) = a₀ ·x⁻ ¹ + a₂ ·x¹ + 0·x² + …

But x⁻ ¹ is not analytic at 0, so this is not a Frobenius series solution around x = 0, and typically
this signals the need for a logarithmic term in the second solution.

We conclude:

One Frobenius solution:

y(x) = a₀ ·x
The second solution involves logarithmic terms (not computed here but theoretically present).

Illustration 2:

Problem:
Use Frobenius method to solve:
x²·y'' + (1 - x)·y' + y = 0

Step 1: Identify p(x), q(x)

Standard form:
x²·y'' + x·p(x)·y' + q(x)·y = 0

We can write:

x²·y'' + (1 - x)·y' + y = 0
Already in suitable form with:

p(x) = (1 - x)/x

x·p(x) = 1 - x

p₀ = 1

q(x) = 1

x²·q(x) = x²

q₀ = 0

Step 2: Indicial equation

r(r - 1) + p₀ ·r + q₀ = 0
r(r - 1) + r = r² - r + r = r² = 0

r=0

Step 3: Assume series solution

y(x) = Σ (n = 0 to ∞) aₙ·xⁿ⁺ ʳ

Use r = 0, so:

y = Σ aₙ·xⁿ
y' = Σ n·aₙ·xⁿ⁻ ¹
y'' = Σ n(n - 1)·aₙ·xⁿ⁻ ²
Step 4: Substitute into ODE
x²·y'' + (1 - x)·y' + y = 0

Substituting:

● x²·y'' = Σ n(n - 1)·aₙ·xⁿ

● (1 - x)·y' = Σ n·aₙ·xⁿ⁻ ¹ - Σ n·aₙ·xⁿ
● y = Σ aₙ·xⁿ

So we get:

Σ n(n - 1)·aₙ·xⁿ + Σ n·aₙ·xⁿ⁻ ¹ - Σ n·aₙ·xⁿ + Σ aₙ·xⁿ = 0

Group all terms:

Σ [n(n - 1)·aₙ - n·aₙ + aₙ]·xⁿ + Σ n·aₙ·xⁿ⁻ ¹ = 0

Shift second sum (change index from n to k = n + 1):

Σ n·aₙ·xⁿ⁻ ¹ = Σ (k = 1 to ∞) k·a_k·x^{k - 1} = Σ (n = 0 to ∞) (n + 1)·a_{n+1}·xⁿ

Now match powers:

Σ [(n(n - 1) - n + 1)·aₙ + (n + 1)·a_{n+1}]·xⁿ = 0

This gives recurrence:

(n² - 2n + 1)·aₙ + (n + 1)·a_{n+1} = 0

Or:

a_{n+1} = -[(n - 1)²·aₙ]/(n + 1)

Step 5: Find coefficients

Let a₀ = C, arbitrary.

● n = 0:
a₁ = -[(0 - 1)²·a₀ ]/1 = -1·C

● n = 1:
a₂ = -[(1 - 1)²·a₁ ]/2 = 0

● n = 2:
a₃ = -[(2 - 1)²·a₂ ]/3 = -1·0/3 = 0

So after a₁ , all terms vanish.

Hence:
y(x) = C·[1 - x]

Illustration 3
Problem:
Solve the differential equation using the Frobenius method:

x²y″ + x(1 − x)y′ − y = 0

Solution:
This is a second-order linear differential equation with a regular singular point at x = 0.

We assume a solution of the form:

y(x) = ∑ (n = 0 to ∞) aₙ xⁿ⁺ ʳ ...(1)

where a₀ ≠ 0, and r is to be determined.

Compute the derivatives:

y′(x) = ∑ (n = 0 to ∞) aₙ (n + r) xⁿ⁺ ʳ⁻ ¹

y″(x) = ∑ (n = 0 to ∞) aₙ (n + r)(n + r − 1) xⁿ⁺ ʳ⁻ ²

Substitute y, y′, and y″ into the original equation:

x²y″ + x(1 − x)y′ − y

= x² ∑ aₙ (n + r)(n + r − 1) xⁿ⁺ ʳ⁻ ²
+ x(1 − x) ∑ aₙ (n + r) xⁿ⁺ ʳ⁻ ¹
− ∑ aₙ xⁿ⁺ ʳ

Simplify each term:

1st term:

x² ∑ aₙ (n + r)(n + r − 1) xⁿ⁺ ʳ⁻ ²
= ∑ aₙ (n + r)(n + r − 1) xⁿ⁺ ʳ

2nd term:

x ∑ aₙ (n + r) xⁿ⁺ ʳ⁻ ¹ − x² ∑ aₙ (n + r) xⁿ⁺ ʳ⁻ ¹
= ∑ aₙ (n + r) xⁿ⁺ ʳ − ∑ aₙ (n + r) xⁿ⁺ ʳ⁺ ¹

3rd term:

− ∑ aₙ xⁿ⁺ ʳ

Now combine all terms:

 ∑ aₙ (n + r)(n + r − 1) xⁿ⁺ ʳ
+ ∑ aₙ (n + r) xⁿ⁺ ʳ
− ∑ aₙ (n + r) xⁿ⁺ ʳ⁺ ¹
− ∑ aₙ xⁿ⁺ ʳ

Group like powers:

∑ xⁿ⁺ ʳ [aₙ ((n + r)(n + r − 1) + (n + r) − 1)]

− ∑ aₙ (n + r) xⁿ⁺ ʳ⁺ ¹

Step 1: Indicial Equation

To find the indicial equation, examine the coefficient of xʳ when n = 0:

a₀ [(r)(r − 1) + r − 1] = 0
a₀ [r² − 1] = 0

Since a₀ ≠ 0, we have:

r² − 1 = 0
r = 1 or r = −1

Step 2: Recurrence Relation

Take the case r = 1 (larger root):

From above:

aₙ ((n + 1)(n) + (n + 1) − 1) = aₙ [(n + 1)(n + 1)]

So the recurrence becomes:

aₙ₊ ₁ = aₙ (n + 1)(n + 1) / ((n + 2)(n + 1))

= aₙ (n + 1) / (n + 2)

This simplifies to:

aₙ₊ ₁ = aₙ (n + 1) / (n + 2)

Step 3: Construct the Series Solution

Let a₀ = arbitrary constant (say, C)

Then:
 a₁ = a₀ (0 + 1)/(1 + 2) = a₀ (1/2)

a₂ = a₁ (1 + 1)/(2 + 2) = (a₀ (1/2)) × (2/4) = a₀ (1/4)

a₃ = a₂ (2 + 1)/(3 + 2) = a₀ (1/4) × (3/5) = a₀ (3/20)

Thus the solution is:

y(x) = x × [a₀ + a₁ x + a₂ x² + a₃ x³ + …]
= x × [C + (C/2)x + (C/4)x² + (3C/20)x³ + …]

Self-Assessment Questions
Instructions:
Use the Frobenius method to find a power series solution about the regular singular point .
Assume the solution of the form
y(x) = ∑ (n = 0 to ∞) aₙ xⁿ⁺ ʳ
and proceed by finding the indicial equation, recurrence relation, and first 3–4 terms of the
series expansion.

Problem 1
x² y″ + 3x y′ − y = 0

Problem 2
x² y″ − x y′ + (x² − 1) y = 0

Problem 3
x² y″ + (2x − x²) y′ − 2y = 0

Problem 4
x² y″ + x y′ + (x² − 4) y = 0
Problem 5
x² y″ − x y′ + y = 0

Problem 6
x² y″ + (4x − x²) y′ − (2 + x) y = 0

Problem 7
x² y″ + (1 − x) y′ − y = 0

Problem 8
x² y″ + (3x + 1) y′ + 2y = 0

Problem 9
x² y″ + x y′ + (x² − 1) y = 0

Problem 10
x² y″ + (2x − 1) y′ − 3y = 0

References
● Boyce, W. E., & DiPrima, R. C. (2017). Elementary Differential Equations and Boundary
Value Problems (11th ed.). Wiley.
● Coddington, E. A. (1955). An Introduction to Ordinary Differential Equations. Prentice
Hall.
● Simmons, G. F. (1972). Differential Equations with Applications and Historical Notes.
McGraw-Hill.
Chapter 4: Legendre Differential Equation –
Theory and Derivation

Learning Objectives
By the end of this chapter, you will be able to:

● Understand the Legendre differential equation and its origin

● Derive the Legendre polynomials and understand their properties
● Apply orthogonality and recurrence relations of Legendre polynomials
● Utilize Legendre functions in physics and engineering applications

4.1. Introduction to Legendre’s Equation

The Legendre differential equation is a second-order linear ordinary differential equation (ODE)
that arises in many physical problems, especially in potential theory, electrostatics, and
spherical harmonics. It naturally emerges when solving Laplace's equation in spherical
coordinates.

The standard form of the Legendre differential equation is:

(1 − x²) d²y/dx² − 2x dy/dx + n(n + 1)y = 0

Where:

● x is the independent variable such that −1 ≤ x ≤ 1

● n is a non-negative integer (n ∈ ℕ ∪ {0})
● y(x) is the dependent function.

4.2. Characteristics
The Legendre equation is:

● A second-order linear ODE.

● Homogeneous.
● With regular singular points at x = ±1.

It is often solved using the Frobenius method or by recognizing its polynomial solution form for
integer n.

4.3. Series Solution

Assume a power series solution of the form:

y(x) = Σₖ ₌₀^∞ aₖ xᵏ

Substitute this series into the Legendre equation and equate powers of x to derive a recurrence
relation for the coefficients aₖ. When n is a non-negative integer, the series truncates, giving
rise to polynomial solutions of degree n. These are called Legendre Polynomials, denoted by
Pₙ(x).

4.4. Legendre Polynomials: Definition

For integer n, the Legendre polynomial Pₙ(x) is defined as the solution to the Legendre
equation that is a polynomial of degree n. These polynomials can also be defined using
Rodrigues’ formula:

Pₙ(x) = (1 / 2ⁿ n!) dⁿ/dxⁿ [ (x² − 1)ⁿ ]

This representation is particularly useful for theoretical derivations and orthogonality proofs.

4.5. Generating Function

The generating function for Legendre polynomials is a compact expression from which all Pₖ(x)
can be derived:

G(x, t) = 1 / √(1 − 2xt + t²) = Σₙ₌ ₀ ^∞ Pₙ(x) tⁿ

Expanding this in a Taylor series in t gives the Legendre polynomials.

4.6. Orthogonality
Legendre polynomials are orthogonal on the interval [−1, 1] with respect to the weight function
w(x) = 1. That is:

∫₋ ₁ ¹ Pₙ(x) Pₙ(x) dx = 0 if n≠m

And for n = m, we have:

∫₋ ₁ ¹ [Pₙ(x)]² dx = 2 / (2n + 1)

This orthogonality property makes Legendre polynomials extremely useful in approximations

and Fourier–Legendre series.
4.7. Recurrence Relations
Legendre polynomials satisfy several recurrence relations. Two important ones are:

● (n + 1) Pₙ₊ ₁ (x) = (2n + 1) x Pₙ(x) − n Pₙ₋ ₁ (x)

● d/dx [Pₙ(x)] = n [Pₙ₋ ₁ (x) − x Pₙ(x)] / (1 − x²)

These relations help in computing higher-order polynomials from lower-order ones.

4.8. Explicit Forms of Legendre Polynomials

Using Rodrigues’ formula:

Pₙ(x) = (1 / 2ⁿ n!) dⁿ/dxⁿ [ (x² − 1)ⁿ ]

The first few Legendre polynomials are:

● P₀ (x) = 1
● P₁ (x) = x
● P₂ (x) = (3x² − 1)/2
● P₃ (x) = (5x³ − 3x)/2
● P₄ (x) = (35x⁴ − 30x² + 3)/8
● P₅ (x) = (63x⁵ − 70x³ + 15x)/8

These are normalized such that Pₙ(1) = 1 for all n

4.9. Parity Property

Legendre polynomials exhibit definite parity:

Pₙ(−x) = (−1)ⁿ Pₙ(x)

● If n is even, Pₙ(x) is even.

● If n is odd, Pₙ(x) is odd.

This is useful when evaluating integrals or solving symmetric boundary value problems.

4.10. Orthogonality Relation

The orthogonality property of Legendre polynomials on [−1, 1] is central to their use in
approximation theory:

∫₋ ₁ ¹ Pₙ(x) Pₙ(x) dx = 0 if n≠m

And:
 ∫₋ ₁ ¹ [Pₙ(x)]² dx = 2 / (2n + 1)

This implies that the set {P₀ (x), P₁ (x), ..., Pₙ(x)} forms an orthogonal basis in the space of
square-integrable functions on [−1, 1].

4.11. Recurrence Relations

Some important recurrence relations are:

(n + 1) Pₙ₊ ₁ (x) = (2n + 1) x Pₙ(x) − n Pₙ₋ ₁ (x)

Useful for computing Pₖ₊ ₁ (x) using Pₖ(x) and Pₖ₋ ₁ (x).

(1 − x²) dPₙ/dx = n x Pₙ(x) − n Pₙ₋ ₁ (x)

Used in evaluating derivatives and solving integrals involving Legendre polynomials.

4.12. Differential Equation Revisited

As a reminder, Legendre polynomials solve the equation:

(1 − x²) d²y/dx² − 2x dy/dx + n(n + 1)y = 0

This ODE has polynomial solutions only when n is a non-negative integer. For non-integer n, the
solutions are Legendre functions of the first and second kind, denoted Pₙ(x) and Qₙ(x),
respectively.

4.13. Application Areas

Legendre polynomials are used in:

● Physics:
○ Potential theory and gravitational fields
○ Electrostatics in spherical coordinates
○ Quantum mechanics (spherical harmonics)

● Numerical Methods
○ Spectral methods for solving PDEs
○ Polynomial approximation and orthogonal expansions

● Engineering:
 ○ Solving Laplace's equation with spherical symmetry
○ Solving boundary value problems

● Mathematics:
○ Fourier–Legendre series expansions
○ Sturm–Liouville theory

Key Terms
● Legendre differential equation
● Legendre polynomials
● Rodrigues’ formula
● Orthogonality
● Generating function
● Recurrence relation
● Parity

Summary / Highlights
● Legendre’s equation appears in problems with spherical symmetry.
● Polynomial solutions exist when n is a non-negative integer.
● Legendre polynomials form an orthogonal basis on [−1, 1].
● Important for physics, numerical methods, and orthogonal expansions.

1. The Legendre differential equation is given by:

A) (1 − x²) y″ − 2x y′ + n(n + 1) y = 0
B) y″ + x y′ + n² y = 0
C) x² y″ + x y′ + (x² − n²) y = 0
D) (1 − x²) y″ + 2x y′ + n(n − 1) y = 0

2. The Legendre polynomial Pₙ(x) is a solution to the Legendre equation for:

A) Real non-integer values of n
B) Integer values of n only
C) Complex values of n
D) Only for n = 0

3. The Rodrigues' formula for Legendre polynomials is:

A) Pₙ(x) = (1 / (2ⁿ n!)) · dⁿ/dxⁿ [(x² − 1)ⁿ]
B) Pₙ(x) = (1 / n!) · dⁿ/dxⁿ [(1 − x²)ⁿ]
C) Pₙ(x) = dⁿ/dxⁿ [(1 + x²)ⁿ]
D) Pₙ(x) = (1 / 2ⁿ) · dⁿ/dxⁿ [(1 − x²)ⁿ]

4. The Legendre polynomials are orthogonal on which interval?

A) [−1, 1]
B) [0, 1]
C) (−∞, ∞)
D) [−π, π]

5. The orthogonality condition of Legendre polynomials is:

A) ∫₋ ¹¹ Pₙ(x) · Pₙ(x) dx = 0, if n = m
B) ∫₋ ¹¹ Pₙ(x) · Pₙ(x) dx = 1
C) ∫₋ ¹¹ Pₙ(x) · Pₙ(x) dx = 0, if n ≠ m
D) ∫₋ ¹¹ Pₙ(x) · Pₙ(x) dx = δₙₙ

6. The degree of the Legendre polynomial Pₙ(x) is:

A) 2n
B) n + 1
C) n
D) n − 1

7. Which of the following is the correct expression for P₀ (x) and P₁ (x)?
A) P₀ (x) = x, P₁ (x) = 1
B) P₀ (x) = 1, P₁ (x) = x
C) P₀ (x) = x², P₁ (x) = x
D) P₀ (x) = 1, P₁ (x) = x²

8. The recurrence relation involving Pₙ(x) is:

A) (n + 1) Pₙ₊ ₁ (x) + (2n + 1) x Pₙ(x) + n Pₙ₋ ₁ (x) = 0
B) (n + 1) Pₙ₊ ₁ (x) − (2n + 1) x Pₙ(x) + n Pₙ₋ ₁ (x) = 0
C) (n + 1) Pₙ₊ ₁ (x) = (2n + 1) x Pₙ(x) − n Pₙ₋ ₁ (x)
D) Pₙ(x) = x Pₙ₋ ₁ (x) − Pₙ₋ ₂ (x)

9. Which of the following Legendre polynomials is correct?

A) P₂ (x) = (1 / 2) (3x² − 1)
B) P₂ (x) = 3x² + 1
C) P₂ (x) = (1 / 2) (x² + 3)
D) P₂ (x) = (1 / 2) (3x − 1)²
10. Which property does not belong to Legendre polynomials?
A) They form a complete orthogonal set on [−1, 1]
B) They are eigenfunctions of a Sturm–Liouville problem
C) They satisfy Bessel’s equation
D) They appear in the solution of Laplace’s equation in spherical coordinates

Answer Key
1. A
2. B
3. A
4. A
5. C
6. C
7. B
8. C
9. A
10. C

Solved Illustrations

Illustration 1

Problem:
Solve the Legendre differential equation for n = 1:

(1 − x²) y″ − 2x y′ + 2y = 0

Solution:
We are given the standard Legendre equation:

(1 − x²) y″ − 2x y′ + n(n + 1) y = 0

Substituting n = 1:

(1 − x²) y″ − 2x y′ + 2 y = 0

Let us try a solution of the form y = P₁ (x), the Legendre polynomial of degree 1.

From known Legendre polynomials:

 P₀ (x) = 1
P₁ (x) = x

Let’s test if y = x satisfies the equation:

We compute the derivatives:

y=x
y′ = 1
y″ = 0

Now substitute into the equation:

(1 − x²)(0) − 2x(1) + 2(x)

= −2x + 2x = 0

Hence, y = x is a solution of the Legendre equation for n = 1.

So the solution is:

y(x) = P₁ (x) = x

Illustration 2:

Problem:
Use Rodrigues’ formula to derive P₂ (x), the second-degree Legendre polynomial.

Rodrigues’ formula:
Pₙ(x) = (1 / (2ⁿ n!)) · dⁿ/dxⁿ [(x² − 1)ⁿ]

Let’s compute for n = 2:

P₂ (x) = (1 / (2² · 2!)) · d²/dx² [(x² − 1)²]

Step 1: Expand the expression inside the derivative:

(x² − 1)² = x⁴ − 2x² + 1

Step 2: Take derivatives:

First derivative:

d/dx [x⁴ − 2x² + 1] = 4x³ − 4x

Second derivative:

d²/dx² [x⁴ − 2x² + 1] = 12x² − 4

Now substitute into Rodrigues’ formula:

P₂ (x) = (1 / (4 · 2)) · (12x² − 4)

= (1 / 8) · (12x² − 4)
= (12x² − 4) / 8
= (3x² − 1) / 2

So, P₂ (x) = (1 / 2)(3x² − 1)

Illustration 3:

Problem:
Solve the Legendre differential equation for n = 2:

(1 − x²) y″ − 2x y′ + 6y = 0

Solution:
We are given:

(1 − x²) y″ − 2x y′ + n(n + 1)y = 0

For n = 2, n(n + 1) = 2 × 3 = 6

Let us try y = P₂ (x), which we know from Rodrigues' formula:

P₂ (x) = (1 / 2)(3x² − 1)

Compute first and second derivatives:

y = (1 / 2)(3x² − 1)
y′ = (1 / 2)(6x) = 3x
y″ = 3

Substitute into the equation:

Left-hand side:

(1 − x²)·3 − 2x·3x + 6·[(1 / 2)(3x² − 1)]

= 3(1 − x²) − 6x² + 3(3x² − 1)

Now compute:

3(1 − x²) = 3 − 3x²

−6x²
3(3x² − 1) = 9x² − 3

Total:
 (3 − 3x²) − 6x² + (9x² − 3)
= (3 − 3x² − 6x² + 9x² − 3)
=0

Hence, y = P₂ (x) = (1 / 2)(3x² − 1) is a solution to the equation.

Illustration 4:

Problem:
Verify that P₁ (x) and P₂ (x) are orthogonal over the interval [−1, 1].

Solution:
Legendre polynomials are orthogonal with respect to the weight function w(x) = 1 on [−1, 1].
This means:

∫ from −1 to 1 of [Pₙ(x)·Pₙ(x)] dx = 0 when m ≠ n

We are given:

P₁ (x) = x
P₂ (x) = (1 / 2)(3x² − 1)

We compute:

I = ∫ from −1 to 1 of [x · (1 / 2)(3x² − 1)] dx

= (1 / 2) ∫ from −1 to 1 of [x(3x² − 1)] dx
= (1 / 2) ∫ from −1 to 1 of [3x³ − x] dx
= (1 / 2) [∫ 3x³ dx − ∫ x dx] from −1 to 1
= (1 / 2) [(3/4)x⁴ − (1/2)x²] from −1 to 1

Evaluate at x = 1:

(3/4)(1)⁴ − (1/2)(1)² = 3/4 − 1/2 = 1/4

Evaluate at x = −1:

(3/4)(−1)⁴ − (1/2)(−1)² = 3/4 − 1/2 = 1/4

So the definite integral is:

(1 / 2)[1/4 − 1/4] = 0

Hence, P₁ (x) and P₂ (x) are orthogonal on [−1, 1].

Self-Assessment Questions

Section A: Basic Application

1. Solve the Legendre differential equation for n = 0:
(1 − x²)y″ − 2x y′ = 0

2. Solve the Legendre differential equation for n = 3:

(1 − x²)y″ − 2x y′ + 12y = 0
Find the polynomial solution.

3. Using Rodrigues' formula, derive P₃ (x).

4. Verify whether the function y(x) = 5x³ − 3x satisfies the Legendre equation for any n.

5. Write the explicit form of P₄ (x) using Rodrigues' formula.

Section B: Orthogonality and Integrals

6. Prove that P₀ (x) and P₁ (x) are orthogonal over the interval [−1, 1].

7. Prove that ∫ from −1 to 1 of [P₂ (x)]² dx = 2 / (2·2 + 1) = 2 / 5

8. Evaluate the integral:

∫ from −1 to 1 of [P₁ (x)·P₃ (x)] dx

9. Let f(x) = x². Expand f(x) in terms of Legendre polynomials up to degree 2.

Section C: Conceptual and Theoretical

10. Show that the general solution to the Legendre equation includes both Pₙ(x) and Qₙ(x),
and explain why Qₙ(x) is usually excluded in physical problems.

11. Derive the recurrence relation:

(n + 1)Pₙ₊ ₁ (x) = (2n + 1)xPₙ(x) − nPₙ₋ ₁ (x)

12. Show that Pₙ(1) = 1 for all n, and explain its significance in boundary value problems.

13. Show that the generating function:

G(x, t) = (1 − 2xt + t²)^(−1/2)
generates Legendre polynomials as the coefficients in the expansion in powers of t.
Prove that the Legendre polynomials satisfy the condition:
d/dx [ (1 − x²) dPₙ(x)/dx ] + n(n + 1)Pₙ(x) = 0

References (APA Style)

● Arfken, G. B., & Weber, H. J. (2005). Mathematical Methods for Physicists (6th ed.).
Academic Press.
● Boyce, W. E., & DiPrima, R. C. (2012). Elementary Differential Equations and Boundary
Value Problems (10th ed.). Wiley.
● Kreyszig, E. (2011). Advanced Engineering Mathematics (10th ed.). Wiley.
● Rainville, E. D. (1971). Special Functions. Macmillan.

End of Chapter 4
Chapter 5: Bessel’s Functions
Learning Objectives
By the end of this chapter, learners will be able to:

● Identify the standard form and origin of Bessel’s equation.

● Apply the Frobenius method to derive Bessel functions.
● Understand the definition and properties of Bessel functions of the first kind.
● Derive and utilize recurrence relations and generating functions.
● Apply orthogonality and asymptotic properties in mathematical and physical problems.

5.1 Introduction to Bessel’s Equation

In applied mathematics and mathematical physics, Bessel’s differential equation arises in
contexts involving cylindrical symmetry such as heat conduction in a circular rod, vibrations of a
circular membrane, and wave propagation in cylindrical structures.

The standard form of Bessel’s equation is:

x² (d²y/dx²) + x (dy/dx) + (x² − n²) y = 0

Where:

● y = y(x) is the unknown function of x,

● n is a real (or integer) constant called the order of the equation.

This is a second-order linear differential equation with variable coefficients.

5.2 Nature and Origin

Bessel’s equation typically results from solving Laplace.equation or the wave equation in
cylindrical coordinates, by using the method of separation of variables. The order n usually
arises from boundary conditions, and corresponds to a mode number or angular dependence in
physical problems

5.3 Frobenius Series Solution

We assume a solution of the form:
y(x) = ∑ (r = 0 to ∞) aᵣ x^(r + s)

Where:

● aᵣ are unknown coefficients to be determined,

● s is an exponent obtained from the indicial equation.

Differentiate term-by-term:

dy/dx = ∑ (r = 0 to ∞) aᵣ (r + s) x^(r + s − 1)

d²y/dx² = ∑ (r = 0 to ∞) aᵣ (r + s)(r + s − 1) x^(r + s − 2)

Substitute into the original Bessel’s equation:

x² (d²y/dx²) + x (dy/dx) + (x² − n²) y = 0

Substitution yields:

x² ∑ aᵣ (r + s)(r + s − 1) x^(r + s − 2)

+ x ∑ aᵣ (r + s) x^(r + s − 1)

+ (x² − n²) ∑ aᵣ x^(r + s)

=0

Now simplify each term:

1. First term becomes:

∑ aᵣ (r + s)(r + s − 1) x^(r + s)

2. Second term becomes:

∑ aᵣ (r + s) x^(r + s)

3. Third term:

∑ aᵣ x^(r + s + 2) − n² ∑ aᵣ x^(r + s)

So the total sum becomes:

∑ aᵣ (r + s)(r + s − 1) x^(r + s)

+ ∑ aᵣ (r + s) x^(r + s)

− ∑ aᵣ n² x^(r + s)

+ ∑ aᵣ x^(r + s + 2) = 0
Combine the first three terms:

∑ aᵣ [(r + s)(r + s − 1) + (r + s) − n²] x^(r + s)

+ ∑ aᵣ x^(r + s + 2) = 0

This gives a recurrence relation for coefficients aᵣ, and the value of s is determined by setting
the coefficient of the lowest power of x to zero — leading to the indicial equation:

s² − n² = 0

s = ±n

Hence, two linearly independent solutions correspond to s = n and s = −n..

5.4 Bessel Function of the First Kind

From the Frobenius method, the Bessel function of the first kind of order n, denoted by
Jₙ(x), is obtained when s = n. Its standard form is:

Jₙ(x) = ∑ (k = 0 to ∞) [ (−1)^k / (k! Γ(k + n + 1)) ] * (x/2)^(2k + n)

Where:

Γ(z) is the Gamma function, with Γ(k + n + 1) = (k + n)! when n is an integer.

Example: J₀ (x)

For n = 0, the series becomes:

J₀ (x) = ∑ (k = 0 to ∞) [ (−1)^k / (k!)² ] * (x² / 4)^k

This function converges for all real values of x and is widely used in practical applications.

5.5 Recurrence Relations

Bessel functions of the first kind, denoted as Jₙ(x), satisfy several important recurrence
relations:

a. First Recurrence Relation

Jₙ₋ ₁ (x) + Jₙ₊ ₁ (x) = (2n / x) × Jₙ(x)

b. Second Recurrence Relation

Jₙ₋ ₁ (x) − Jₙ₊ ₁ (x) = 2 × dJₙ(x)/dx

This can also be rearranged to write the derivative explicitly:

c. Derivative Formula

dJₙ(x)/dx = (1/2) × [Jₙ₋ ₁ (x) − Jₙ₊ ₁ (x)]

These recurrence relations are useful for simplifying expressions and solving differential
equations involving Bessel functions.

and solving recurrence-based PDEs.

5.6 Orthogonality Property

If αₙ,ₙ and αₙ,ₙ are the m-th and k-th zeros of Jₙ(x) respectively, and R is a positive
constant, then the orthogonality condition is:

∫₀ᴿ x × Jₙ((αₙ,ₘ × x) / R) × Jₙ((αₙ,ₖ × x) / R) dx = 0 if m ≠ k

This property is important in expanding functions using the Bessel series in cylindrical problems.

5.7 Generating Function

The generating function for Bessel functions of integer order is given by:

G(x, t) = exp[(x / 2) × (t − 1/t)]

The power series expansion of G(x, t) is:

G(x, t) = Σ (tⁿ × Jₙ(x)) where the sum is from n = −∞ to ∞

This generating function is a compact expression that produces all Bessel functions as
coefficients of powers of t.

5.8 Integral Representation

The Bessel function of the first kind can also be represented using an integral formula:

Jₙ(x) = (1 / π) × ∫₀ ^π cos(nθ − x × sinθ) dθ

This is valid for real x and integer n, and is especially useful in analytical derivations and
numerical calculations.

5.9 Asymptotic Behavior

As x → ∞, the Bessel function Jₙ(x) behaves like a decaying cosine wave:

Jₙ(x) ≈ √(2 / (πx)) × cos(x − (nπ / 2) − (π / 4))

This approximation is helpful in studying the behavior of solutions at large distances or time
scales.

5.10 Zeros of Bessel Functions

The Bessel function Jₙ(x) has infinitely many real and positive zeros.

Let αₙ,ₙ denote the m-th positive zero of Jₙ(x).

These zeros are crucial in boundary value problems in cylindrical domains, especially in solving
PDEs with radial symmetry.

Summary
● Bessel’s equation arises in cylindrical systems.
● Solutions are expressed via Frobenius method.
● is the Bessel function of the first kind.
● Bessel functions have useful orthogonality and recurrence properties.
● They are essential in boundary value problems with circular symmetry.

Check Your Understanding

1. The Bessel differential equation is:

A) x²y'' + xy' + (x² − n²)y = 0

B) x²y'' − xy' + (x² − n²)y = 0
C) xy'' + x²y' + (n² − x²)y = 0
D) y'' + (1/x)y' + y = 0
2. The Bessel function of the first kind of order n, denoted by Jₙ(x), is defined by the series:

A) Jₙ(x) = ∑ₙ₌ ₀ ^∞ [ (−1)^k * (x²)^k ] / [k! * Γ(n + k + 1)]

B) Jₙ(x) = ∑ₙ₌ ₀ ^∞ [ (−1)^k * (x/2)^(2k+n) ] / [k! * Γ(n + k + 1)]
C) Jₙ(x) = ∑ₙ₌ ₀ ^∞ [ x^k ] / [k! * Γ(n − k + 1)]
D) Jₙ(x) = ∑ₙ₌ ₀ ^∞ [ x^(n+2k) ] / [Γ(k+1)]

3. Which of the following is a recurrence relation for Bessel functions?

A) Jₙ₊ ₁ (x) + Jₙ₋ ₁ (x) = (2n/x) Jₙ(x)

B) Jₙ₊ ₁ (x) − Jₙ₋ ₁ (x) = (2/x) Jₙ(x)
C) xJₙ₋ ₁ (x) + xJₙ₊ ₁ (x) = 2n Jₙ(x)
D) Jₙ(x) = (x/2)^n / Γ(n+1)

4. The orthogonality relation for Bessel functions is given by:

A) ∫₀ ^∞ Jₙ(x) dx = 0
B) ∫₀ ^∞ x Jₙ(αₙx) Jₙ(αₙx) dx = 0 for m ≠ n
C) ∫₀ ^∞ x² Jₙ(x) dx = δₙₙ
D) ∫₀ ^∞ Jₙ(x) Jₙ(x) dx = π/2

5. The Bessel function of the second kind is denoted by:

A) Yₙ(x)
B) Jₙ(x)
C) Iₙ(x)
D) Kₙ(x)

6. The Bessel function of order zero, J₀ (x), satisfies:

A) J₀ ''(x) + J₀ (x) = 0
B) x²J₀ ''(x) + xJ₀ '(x) + x²J₀ (x) = 0
C) x²J₀ ''(x) + xJ₀ '(x) + xJ₀ (x) = 0
D) x²J₀ ''(x) − xJ₀ '(x) + x²J₀ (x) = 0
7. The generating function for Bessel functions of the first kind is:

A) e^(xt) = ∑ₙ₌ ₀ ^∞ Jₙ(x)tⁿ

B) e^(x(t − 1/t)/2) = ∑ₙ₌ ₋ ∞^∞ Jₙ(x) tⁿ
C) sin(xt) = ∑ₙ₌ ₀ ^∞ Jₙ(x)tⁿ
D) e^(xt²) = ∑ₙ₌ ₀ ^∞ Jₙ(x)tⁿ

8. For small values of x, the Bessel function Jₙ(x) behaves approximately as:

A) Jₙ(x) ≈ xⁿ
B) Jₙ(x) ≈ sin(x)
C) Jₙ(x) ≈ x^(−n)
D) Jₙ(x) ≈ 1/x

Answer Key
1. A
2. B
3. A
4. B
5. A
6. B
7. B
8. A

Solved Illustration 1:
Problem:
Show that the Bessel function of the first kind of order zero,
J₀ (x) = ∑ₙ₌ ₀ ^∞ [ (−1)^k * (x²/4)^k ] / (k!)²
satisfies the Bessel differential equation:
x²y'' + x y' + x²y = 0

Solution:
We are given:
y = J₀ (x) = ∑ₙ₌ ₀ ^∞ [ (−1)^k * (x²/4)^k ] / (k!)²

Let us compute the first and second derivatives:

First derivative:
y' = d/dx [J₀ (x)]
= ∑ₙ₌ ₁ ^∞ [ (−1)^k * k * 2x * (x²/4)^(k−1) ] / (k!)²
= ∑ₙ₌ ₁ ^∞ [ (−1)^k * x * (x²/4)^(k−1) ] / [ (k−1)! * k! ]

Second derivative:
y'' = d²/dx² [J₀ (x)]
= ∑ₙ₌ ₁ ^∞ [ (−1)^k * ((2k−1)x² + 2(k−1)) * (x²/4)^(k−2) ] / (k!)²
(Derivative gets complicated, but we proceed symbolically)

Now substitute into the equation:

x²y'' + x y' + x²y = 0

If we do the substitution term by term using power series (as in standard Bessel analysis), the
result indeed simplifies to 0, thereby verifying that J₀ (x) satisfies the differential equation.

Conclusion:
Hence, J₀ (x) is a solution of the Bessel equation of order zero.

Solved Illustration 2:
Problem:
Find the first four non-zero terms of the Bessel function of the first kind of order one:
J₁ (x) = ∑ₙ₌ ₀ ^∞ [ (−1)^k * (x/2)^(2k+1) ] / [k! * Γ(k + 2)]

Solution:
Use the definition:
J₁ (x) = ∑ₙ₌ ₀ ^∞ [ (−1)^k * (x/2)^(2k+1) ] / [k! * (k + 1)!]

Let us compute terms up to k = 3:

● For k = 0:
Term₀ = (−1)^0 * (x/2)^1 / (0! * 1!) = x/2

● For k = 1:
Term₁ = (−1)^1 * (x/2)^3 / (1! * 2!) = − x³ / (2³ * 2) = − x³ / 16

● For k = 2:
Term₂ = (−1)^2 * (x/2)^5 / (2! * 3!) = x⁵ / (2⁵ * 12) = x⁵ / 384

● For k = 3:
Term₃ = (−1)^3 * (x/2)^7 / (3! * 4!) = − x⁷ / (128 * 144) = − x⁷ / 18432
Thus,
J₁ (x) ≈ x/2 − x³/16 + x⁵/384 − x⁷/18432

Conclusion:
The series expansion of J₁ (x) up to four terms is:
J₁ (x) ≈ x/2 − x³/16 + x⁵/384 − x⁷/18432

Solved Illustration 3:
Problem:
Given that:
J₀ (1) ≈ 0.7652 and J₁ (1) ≈ 0.4401,
verify the recurrence relation:
J₀ (x) = (2/x) * J₁ (x) − J₂ (x),
and find an approximate value of J₂ (1).

Solution:
We are given:
J₀ (1) = 0.7652
J₁ (1) = 0.4401
x=1

Rewriting the recurrence:

J₂ (1) = (2/1) * J₁ (1) − J₀ (1)
= 2 * 0.4401 − 0.7652
= 0.8802 − 0.7652
= 0.1150

Conclusion:
The approximate value of J₂ (1) is 0.1150, which satisfies the recurrence relation.

Solved Illustration 4:
Problem:
Show that Bessel functions of the first kind of order zero are orthogonal with respect to the
weight function xx, i.e.,
∫₀ ^a x J₀ (αₙx/a) J₀ (αₙx/a) dx = 0
for m ≠ n, where αₙ and αₙ are the m-th and n-th positive roots of J₀ (x) = 0.
Solution:
Let αₙ and αₙ be distinct positive zeros of J₀ (x), i.e.,
J₀ (αₙ) = 0 and J₀ (αₙ) = 0.

The orthogonality property of Bessel functions states:

∫₀ ^a x J₀ (αₙx/a) J₀ (αₙx/a) dx = 0 when m ≠ n.

This is a standard result of Sturm-Liouville theory applied to Bessel's differential equation,

where J₀ (x) is an eigenfunction and the zeros αₙ are eigenvalues.

We do not compute it explicitly here due to complexity, but the result is proven in advanced
texts using integration by parts and properties of Bessel's differential operator.

Conclusion:
The Bessel functions J₀ (αₙx/a) are orthogonal over [0, a] with weight x when m ≠ n.

Solved Illustration 5:
Problem:
Evaluate the integral:
∫₀ ^a x J₁ (kx) dx

Solution:
We use the known result:
∫₀ ^a x J₁ (kx) dx = a² J₂ (ka) / (2k)

Let us take an example:

Let a = 1, and k = 2. Then:

∫₀ ^1 x J₁ (2x) dx = (1)² * J₂ (2) / (2 * 2)

= J₂ (2) / 4

From tables or software:

J₂ (2) ≈ 0.3528

Thus:
∫₀ ^1 x J₁ (2x) dx ≈ 0.3528 / 4 = 0.0882

Conclusion:
∫₀ ^1 x J₁ (2x) dx ≈ 0.0882
 Unsolved Problem Set: Bessel Functions
Problem 1:
Show that the Bessel function Jn(x)J_n(x) satisfies the identity:
x Jₙ₋ ₁ (x) + x Jₙ₊ ₁ (x) = 2n Jₙ(x)

Problem 2:
Derive the recurrence relation:
Jₙ₋ ₁ (x) - Jₙ₊ ₁ (x) = 2 Jₙ′(x)
where Jₙ′(x) denotes the derivative of Jₙ(x) with respect to x.

Problem 3:
Prove that the function
y(x) = xⁿ
is not a solution to Bessel’s differential equation for any non-zero order n.

Problem 4:
Use the method of Frobenius to derive the series solution of Bessel’s differential equation of
order zero:
x²y″ + x y′ + x² y = 0

Problem 5:
If
y(x) = J₁ (x),
show that
d/dx [x J₁ (x)] = x J₀ (x)

Problem 6:
Evaluate the integral:
∫₀ ^1 x J₀ (kx) dx
where k is a constant.

Problem 7:
Let α₁ , α₂ , ... be the positive zeros of J₀ (x). Show that the functions
ϕₙ(x) = J₀ (αₙ x)
are orthogonal over the interval [0, 1] with respect to the weight function x.

Problem 8:
Verify that
Jₙ(−x) = (−1)ⁿ Jₙ(x)
and hence determine the parity (evenness or oddness) of Bessel functions.

Problem 9:
Find the value of the Wronskian
 W[Jₙ(x), Yₙ(x)]
where Jₙ(x) and Yₙ(x) are the Bessel functions of the first and second kind respectively.

Problem 10:
Using recurrence relations, find J₄ (x) in terms of J₃ (x) and J₂ (x).

References
● Boyce, W. E., & DiPrima, R. C. (2017). Elementary Differential Equations and Boundary
Value Problems. Wiley.
● Zill, D. G. (2018). A First Course in Differential Equations with Modeling Applications.
Cengage Learning.
● Arfken, G. B., Weber, H. J., & Harris, F. E. (2012). Mathematical Methods for Physicists.
Academic Press.
Module 4: Boundary Value Problems
Chapter 1: Introduction to Boundary Value Problems

Learning Objectives

By the end of this chapter, the learner will be able to:

● Define and distinguish boundary value problems (BVPs) from initial value problems
(IVPs).
● Understand the structure of second-order linear BVPs and the concept of well-
posedness.
● Classify BVPs based on regularity, type of equation, boundary conditions, and linearity.
● Apply analytical and numerical methods to solve BVPs.
● Recognize applications of BVPs in physical sciences and engineering.
● Interpret the role of eigenvalue problems and Sturm–Liouville theory in the context of
BVPs.

4.1 Introduction to Boundary Value Problems (BVPs)

Boundary Value Problems (BVPs) are a critical class of problems in the study of differential
equations, especially in mathematical physics and engineering, where a differential equation is
solved subject to conditions (constraints) specified at multiple points, frequently at the domain
boundary.

BVPs specify the solution behaviour at the boundaries of the domain, usually at two different
values of the independent variable, as opposed to Initial Value Problems (IVPs), where all the
information about the solution is provided at a single point (often the initial point).

Definition:

A Boundary Value Problem (BVP) is a differential equation together with a set of additional
constraints, called boundary conditions, that the solution must satisfy at the endpoints of the
interval.

A general second-order linear boundary value problem is of the form:

(1)
d²y/dx² + p(x) dy/dx + q(x) y = f(x),
for a ≤ x ≤ b,
subject to boundary conditions such as:
y(a) = α, y(b) = β

Here:

● p(x), q(x), f(x) are known functions,

● [a, b] is the domain (interval),

● α, β are known constants (boundary values).

Key Concepts in BVPs:

1. Well-posedness:
A BVP is considered well-posed if:
○ A solution exists,
○ The solution is unique,
○ The solution depends continuously on the boundary data.

2. Types of Boundary Conditions:

BVPs can have various types of boundary conditions:
○ Dirichlet boundary conditions (values of the function at the endpoints):
y(a) = A, y(b) = B
○ Neumann boundary conditions (values of the derivative at the endpoints):
dy/dx(a) = A, dy/dx(b) = B
○ Robin (mixed) boundary conditions (linear combinations of values and
derivatives):
a₁y(a) + b₁dy/dx(a) = A
a₂y(b) + b₂dy/dx(b) = B

3. Linear vs. Nonlinear BVPs:

○ If the differential operator and boundary conditions are linear, the BVP is linear.
○ If nonlinearity exists in either the differential equation or the boundary conditions,
it becomes a nonlinear BVP.

4.2 Difference Between Initial and Boundary Value Problems

Feature Initial Value Problem (IVP) Boundary Value Problem
(BVP)
 Specification of data At a single point (x = a) At multiple points (x = a, x =
b)

Solution method Integration, IVP solvers Eigenfunction expansion,

finite difference

Applications Population, circuits, Heat conduction, vibrations,

mechanics quantum models

Type of solution Forward (causal) Global (domain-constrained)

4.3 Applications of Boundary Value Problems

Boundary value problems arise naturally in many physical systems governed by differential
equations. Some examples include:

● Heat conduction in a rod:

Described by the second-order partial differential equation, reduced to a BVP in steady
state.
Example:
d²T/dx² = 0, with T(0) = T₀, T(L) = T_L

● Vibrating string or membrane:

Wave equation leads to BVPs when looking for normal modes or steady-state solutions.

● Electrostatics and Fluid Dynamics:

Laplace's or Poisson's equations with potential or pressure boundary values.

● Quantum mechanics:
Schrödinger’s equation as a BVP for the energy eigenvalue problem.

● Structural analysis:
Bending of beams and plates modeled via differential equations with deflection and
slope at ends.

4.4 Graphical Interpretation

Finding a function that "fits" the shape given by the differential equation while also meeting two
or more "pins" or constraints at the boundary is how a BVP is represented. In IVPs, on the
other hand, the function merely "grows" from a known beginning shape.
4.5 Classification of BVPs

Boundary value problems can be classified in various ways based on their structure and
properties.

(a) Regular vs Singular BVPs

1. Regular BVP:
A boundary value problem is called regular if:
○ The functions p(x), q(x), and f(x) are continuous over the closed interval
[a, b].
○ The coefficient of the highest-order derivative (usually 1 or non-zero function)
does not vanish within [a, b].
○ The boundary conditions are finite and well-defined at x = a and x = b.
Example:
d²y/dx² + y = sin(x), with y(0) = 0, y(π) = 0

2. Singular BVP:
A BVP is called singular if:
○ One or more of the functions p(x), q(x), or the coefficient of the highest
derivative is discontinuous or unbounded on [a, b].
○ The domain includes points where the equation is not defined (such as division
by zero or infinite behavior).
Example:
d/dx(x dy/dx) + y = 0, with y(0) and y(1) given.
This is singular at x = 0, where x dy/dx becomes undefined.

4.6 Methods of Solving BVPs

There are several analytical and numerical methods used to solve BVPs. The method is chosen
depending on their complexity and form.

(a) Analytical Methods

1. Separation of Variables:
This method is used when the problem can be decomposed into independent functions
of x and t (or other variables). Often leads to eigenvalue problems.

2. Method of Eigenfunction Expansion:

This method expands the solution as a series of orthogonal eigenfunctions derived from
 a related Sturm–Liouville problem.

3. Green’s Function Method:

Constructs a special function called Green’s function to express the solution in integral
form.

4. Transform Methods:
Techniques like the Laplace or Fourier transform convert the differential equation into an
algebraic one, which is easier to solve.

(b) Numerical Methods

When exact solutions are difficult or impossible to find, numerical techniques are used:

1. Shooting Method:
○ Converts the BVP into an initial value problem (IVP).
○ Guess initial slope and integrate.
○ Adjust the guess iteratively until the boundary condition at the other endpoint is
satisfied.

Idea:

Solve d²y/dx² = f(x, y, y′)

Guess y′(a) = s, solve IVP: y(a) = α, y′(a) = s
Then adjust s such that y(b) = β

2. Finite Difference Method (FDM):

○ Replace derivatives with difference quotients.
○ Leads to a system of linear equations.
○ Solves for approximated values of the function at grid points.

Example:
Approximate d²y/dx² using:
(y_{i+1} - 2y_i + y_{i-1}) / h² ≈ d²y/dx²
Build a system and solve.

3. Finite Element Method (FEM):

○ Divide the domain into small elements.
○ Construct local approximations (piecewise polynomials).
○ Assemble a global system using variational methods.
○ Common in engineering applications for irregular domains.
4.7 Eigenvalue Problems in BVPs

Many BVPs arise as eigenvalue problems, especially in physics where boundary conditions
lead to quantized solutions (e.g., vibrations, quantum mechanics).

Definition:

An eigenvalue problem is a BVP of the form:

L[y] = λ y,
with boundary conditions at x = a and x = b.

Here:

● L is a differential operator (often self-adjoint),

● λ is an unknown constant (eigenvalue),
● Nontrivial solutions y(x) that satisfy both the differential equation and boundary
conditions are called eigenfunctions.

Example:

d²y/dx² + λ y = 0,
with y(0) = 0, y(π) = 0

The general solution:

y(x) = A sin(√λ x) + B cos(√λ x)

Apply boundary conditions:

● y(0) = 0 ⇒ B = 0
● y(π) = 0 ⇒A sin(√λ π) = 0
⇒ √λ π = nπ, so λ = n², n ∈ ℕ

Thus:

● Eigenvalues: λₙ = n²
● Eigenfunctions: yₙ(x) = sin(n x)

4.8 Properties of Eigenvalue Problems

 ● Orthogonality:
If yₙ and yₘ are eigenfunctions corresponding to distinct eigenvalues λₙ ≠ λₘ, then:
∫ₐᵇ yₙ(x) yₘ(x) dx = 0

● Completeness:
The set of eigenfunctions forms a complete basis. Any function which satisfies the
boundary conditions can be expanded as a linear combination of eigenfunctions.

● Sturm–Liouville Form:
Many eigenvalue problems can be written as:
d/dx[p(x) dy/dx] + [λ w(x) - q(x)] y = 0
with boundary conditions.
This is the general form of a Sturm–Liouville problem, central to understanding BVPs

4.9 Further Classification of BVPs

(A) By Type of Differential Equation

Ordinary Differential Equations (ODEs):
These BVPs are defined on a single independent variable (typically time or one-dimensional
space). For example, consider the equation:
d²y/dx² + λy = 0

with boundary conditions:

y(0) = 0, y(π) = 0

This describes many physical systems like a vibrating string fixed at both ends, and is an ODE-
based BVP.

Partial Differential Equations (PDEs):

These problems involve two or more independent variables (typically space and time) and are
common in multidimensional systems such as heat conduction, wave motion, or fluid flow. For
example:

∂²u/∂x² + ∂²u/∂y² = 0

it is the two-dimensional Laplace’s equation with boundary values specified on the boundary of
a domain.
(B) By Boundary Conditions

Boundary conditions define how the solution behaves at the edges of the domain. The three
classical boundary conditions are:

Dirichlet Conditions (First Kind):

The values of the function are fixed at the boundaries.

Example:
y(0) = 0, y(1) = 2

This often represents physically clamped or fixed systems—e.g., the ends of a string nailed
down.

Neumann Conditions (Second Kind):

The values of the derivative (typically spatial gradient) are fixed.

Example:
dy/dx (0) = 0, dy/dx (1) = 1

This models situations like insulated boundaries in heat conduction, where the flux is constant
or zero.

Robin (Mixed) Conditions (Third Kind):

A linear combination of function and derivative values is specified.

Example:
a·y(0) + b·dy/dx (0) = c

These arise in convection problems and other physical systems where there is interaction
between a system and its surroundings.

(C) By Linearity

Linear BVPs:
Both the differential equation and the boundary conditions are linear. These are easier to
analyze and solve using superposition, eigenfunction expansions, and other analytical tools.

Nonlinear BVPs:
These involve nonlinear operators or boundary conditions. For instance:
d²y/dx² + y³ = 0, y(0) = 0, y(1) = 1

Such equations require special numerical methods and are more complex due to issues of non-
uniqueness, bifurcations, and sensitivity.

4.10 BVP vs IVP: Conceptual Contrast

Aspect IVP BVP

Data Location Single point Multiple boundary points

Nature of Solution Forward evolving Spatially constrained

Existence/Uniqueness Often guaranteed Depends on well-posedness

Methods Runge–Kutta, integration Sturm–Liouville, Green's,

FEM/FDM

Examples Growth models, motion Heat flow, fluid pressure,

quantum

Insight:
IVPs describe evolution over time; BVPs define configuration over space.

4.11 Significance and Real-World Applications

(i) Heat Transfer

In a rod of length LL, with temperatures fixed at both ends, the steady-state heat equation
becomes:

d²T/dx² = 0, T(0) = T₀ , T(L) = T₁

This BVP determines how temperature varies spatially in a non-evolving (steady) system.

(ii) Vibrations and Acoustics

The vibration of a stretched string of length L fixed at both ends is governed by:

d²y/dx² + (λ/c²)·y = 0, y(0) = 0, y(L) = 0

Solving this BVP gives the natural frequencies (eigenvalues) and vibration modes
(eigenfunctions), which are crucial in designing musical instruments, bridges, and buildings.

(iii) Electrostatics

Laplace’s equation ∇ 2ϕ=0∇ ²ϕ = 0, subject to boundary potential values, defines how electric
potential behaves inside a closed domain. This is vital in electric field modeling, capacitor
design, and electromagnetic shielding.

(iv) Quantum Mechanics

The Schrödinger equation is a second-order differential equation where the boundary conditions
reflect physical confinement (e.g., a particle in a box):

d²ψ/dx² + (2m/ħ²)(E − V(x))ψ = 0

Here, solving the BVP gives the energy levels of the particle and wavefunctions ψ(x)ψ(x).

(v) Fluid Flow and Mechanics

Incompressible fluid velocity profiles (e.g., Poiseuille flow in a pipe) arise from solving BVPs
using the Navier–Stokes equations with boundary conditions like no-slip walls.

4.12 Check Your Progress

1. A boundary value problem (BVP) differs from an initial value problem (IVP) because:
A) BVP requires only initial data for the solution
B) BVP specifies the solution at more than one point
C) BVP applies only to nonlinear equations
D) BVP cannot be solved numerically

2. Consider the differential equation:

d²y/dx² + λy = 0, with boundary conditions y(0) = 0 and y(L) = 0.
This problem is an example of:
A) Initial value problem
B) Dirichlet boundary value problem
C) Neumann boundary value problem
D) Robin boundary value problem

3. In a Sturm–Liouville problem, the differential equation takes the form:

d/dx[p(x) dy/dx] + [λw(x) − q(x)] y = 0
Here, λ is known as:
A) A constant
B) A forcing term
C) An eigenvalue
D) A Lagrange multiplier

4. The Green’s function G(x, ξ) for a linear differential operator L satisfies:

A) L[G(x, ξ)] = 0
B) L[G(x, ξ)] = δ(x − ξ)
C) G(x, ξ) = δ′(x − ξ)
D) L[G(x, ξ)] = 1

5. For the Sturm–Liouville problem, the eigenfunctions corresponding to distinct eigenvalues

are:
A) Linearly dependent
B) Orthogonal with respect to the weight function w(x)
C) Identical
D) Nonexistent

6. Which of the following is not true for a regular Sturm–Liouville problem?

A) All eigenvalues are real
B) Eigenfunctions corresponding to distinct eigenvalues are orthogonal
C) Eigenfunctions form a complete basis
D) Eigenvalues can be complex in general

7. In the boundary value problem

d²y/dx² + y = 0, with y(0) = 0 and y(π) = 0,
what are the eigenvalues λ?
A) λ = n², where n = 1, 2, 3,...
B) λ = nπ, where n = 1, 2, 3,...
C) λ = n²π², where n = 1, 2, 3,...
D) λ = −n²π², where n = 1, 2, 3,...

8. The application of Green’s function is primarily used to:

A) Solve algebraic equations
B) Solve non-homogeneous linear differential equations
C) Solve nonlinear boundary value problems
D) Solve optimization problems

9. The function δ(x − ξ) used in the definition of Green’s function is:

A) A continuous function
B) A step function
C) A Dirac delta function
D) A periodic function

10. For the boundary condition

a₁ y(a) + a₂ y′(a) = 0,
this is an example of:
A) Dirichlet condition
B) Neumann condition
C) Robin (mixed) boundary condition
D) Periodic boundary condition

Answer Key

1. B
2. B
3. C
4. B
5. B
6. D
7. C
8. B
9. C
10. C
4.13 Summary

Boundary Value Problems are foundational in modeling equilibrium and steady-state behavior
across disciplines. Through classification, analytical and numerical solution techniques, and the
role of eigenvalue problems, BVPs emerge as indispensable tools in scientific computation,
engineering design, and mathematical physics. Sturm–Liouville theory and Green’s functions
provide theoretical strength, while FDM and FEM allow computational implementation.

4.14 Key Terms

● Boundary Value Problem (BVP): A differential equation with constraints at more than
one point.
● Well-posedness: A problem that is solvable, has a unique solution, and is stable under
small changes.
● Dirichlet Condition: Fixed value of the function at boundary.
● Neumann Condition: Fixed derivative at boundary.
● Robin Condition: Linear combination of function and derivative at boundary.
● Sturm–Liouville Problem: A second-order linear BVP with spectral parameter λ.
● Eigenvalue: A scalar λ for which a non-trivial solution exists.
● Green’s Function: A special function used to solve nonhomogeneous BVPs.
● Orthogonality: Mutual independence of eigenfunctions under inner product.
● Shooting Method: Converts BVP to IVP iteratively.

4.15 Self-Assessment Questions

ₙ Short Answer

1. Define a Boundary Value Problem.

2. What are the key differences between BVP and IVP?
3. Give one example of a physical system modeled by a BVP.
4. Explain what makes a BVP regular.
5. What is a Sturm–Liouville problem?
6. Write the finite difference approximation for d²y/dx².
7. What are Dirichlet boundary conditions?
8. Why are eigenfunctions important in solving PDEs?

ₙ Long Answer

1. Compare and contrast the classification of BVPs based on linearity and boundary
conditions.
2. Explain the method of Green’s functions for solving linear BVPs.
 3. Derive eigenvalues and eigenfunctions for the problem d²y/dx² + λy = 0, with y(0) = 0,
y(π) = 0.
4. Discuss the applications of BVPs in quantum mechanics and fluid flow.
5. Describe the shooting method and finite difference method for solving BVPs.
6. Explain Sturm–Liouville theory and its significance in solving BVPs.
7. Describe the real-world implications of well-posedness in boundary value problems.

4.16 Case Study: Heat Conduction in a Metal Rod

Background:
A rod of length L is insulated on the sides and has fixed temperatures T₀ and T₁ at x = 0 and x
= L respectively. The temperature distribution is governed by:

d²T/dx² = 0, T(0) = T₀ , T(L) = T₁

Solution:
The general solution of d²T/dx² = 0 is:

T(x) = Ax + B

Applying boundary conditions:

● T(0) = B = T₀
● T(L) = AL + B = T₁ ⇒ A = (T₁ − T₀ )/L

So:

T(x) = T₀ + (T₁ − T₀ )x/L

Significance:
This solution represents a linear temperature profile — a classical example of a regular, linear
BVP.

Assessment Questions:

1. What type of boundary conditions are used in this problem?

2. Why is the solution linear?
3. What would change if the rod lost heat to the environment?

4.17 References
● Boyce, W. E., & DiPrima, R. C. (2017). Elementary Differential Equations and Boundary
Value Problems (11th ed.). Wiley.
● Haberman, R. (2012). Applied Partial Differential Equations (5th ed.). Pearson.
● Zill, D. G. (2012). A First Course in Differential Equations with Modeling Applications
(10th ed.). Cengage Learning.
● Simmons, G. F. (1972). Differential Equations with Applications and Historical Notes.
McGraw-Hill.
Chapter 2: Sturm–Liouville Problems
2.1 Learning Objectives
● Understand the historical and mathematical context of Sturm–Liouville theory.
● Recognize the general form and structure of Sturm–Liouville differential equations.
● Explain self-adjointness and the role of inner products.
● Identify and apply appropriate boundary conditions.
● Analyze eigenvalue problems and their properties.
● Explore real-world applications of Sturm–Liouville problems.

2.2 Historical and Mathematical Context

The 19th-century French mathematicians Charles Sturm and Joseph Liouville laid the
foundation for the Sturm–Liouville theory. This theory provides a fundamental framework for the
analysis of linear second-order differential equations with variable coefficients and is essential in
solving boundary value problems (BVPs) in mathematical physics. It unifies various special
functions (e.g., Legendre, Bessel, Hermite) and introduces powerful tools for spectral analysis,
orthogonal expansions, and separation of variables in partial differential equations (PDEs).

2.3 General Form of Sturm–Liouville Differential Equation

A second-order linear differential equation of the Sturm–Liouville form is given by:

d/dx [ p(x) * dy/dx ] + [ λ * w(x) - q(x) ] * y = 0

Or, equivalently written:

( d/dx )[ p(x) * dy/dx ] + [ λ * w(x) - q(x) ] * y(x) = 0

Where:

● y(x) is the unknown function,

● λ is the eigenvalue parameter,
● p(x), q(x), and w(x) are known real-valued functions defined on a closed interval [a,
b],
● p(x) > 0, w(x) > 0 on [a, b] for physical relevance.

This is often written compactly as:

L[y] = λ * w(x) * y

Where L is the Sturm–Liouville operator:

L[y] = - d/dx [ p(x) * dy/dx ] + q(x) * y

2.4 Self-Adjoint Form

A key property of Sturm–Liouville operators is that they are self-adjoint (or Hermitian) under an
appropriate inner product. That is:

⟨ L[y], z⟩ = ⟨ y, L[z]⟩

with respect to the inner product defined by:

⟨f, g⟩ = ∫ₐ ᵇ f(x) * g(x) * w(x) dx

This property ensures that the eigenvalues λ are real and the eigenfunctions corresponding to
distinct eigenvalues are orthogonal with respect to the weight function w(x).

2.5 Boundary Conditions Boundary

● To fully specify a Sturm–Liouville problem, the differential equation is paired with
homogeneous boundary conditions at the endpoints x = a and x = b. These may
take several forms:
● Common types of boundary conditions:
● Dirichlet (fixed ends):

y(a) = 0, y(b) = 0
● Neumann (free ends):

y'(a) = 0, y'(b) = 0
● Robin (mixed):

α₁ * y(a) + β₁ * y'(a) = 0
● α₂ * y(b) + β₂ * y'(b) = 0
● These conditions ensure that the operator L remains self-adjoint and that the
eigenfunctions form an orthogonal basis under the inner product.
2.6 Eigenvalue Problems
A Sturm–Liouville problem is thus an eigenvalue problem where we seek values of λ (called
eigenvalues) such that the differential equation has non-trivial solutions y(x) (called
eigenfunctions) satisfying the boundary conditions.

The collection of all such eigenvalues forms the spectrum of the problem. Each eigenvalue λₙ
has a corresponding eigenfunction yₙ(x), and under suitable conditions:

{ y₁ (x), y₂ (x), y₃ (x), ... }

form an orthogonal basis in the Hilbert space L²[a, b] with weight w(x)

2.7 Key Theoretical Properties

The Sturm–Liouville differential equation is of the form:

L[y] = d/dx (p(x) · dy/dx) + q(x) · y = -λ · w(x) · y

with prescribed boundary conditions on the interval [a, b]. The properties of the solutions
(eigenfunctions) and associated values (eigenvalues) are foundational to many
applications in mathematical physics and engineering. Below is an elaboration of each core
property:

Property 1: Orthogonality

Description: Eigenfunctions corresponding to different eigenvalues are orthogonal with respect

to the weight function w(x).

Let yₘ(x) and yₙ(x) be two eigenfunctions of the Sturm–Liouville problem with distinct
eigenvalues λₘ ≠ λₙ . Then the orthogonality relation is:

∫ₐᵇ w(x) · yₘ(x) · yₙ(x) dx = 0

This means that the eigenfunctions are mutually perpendicular under a weighted inner product
space. Orthogonality is a consequence of the self-adjointness of the differential operator L.

Property 2: Real Eigenvalues

Description: All eigenvalues λ are guaranteed to be real numbers due to the self-adjoint
nature of the operator L.

If the Sturm–Liouville operator is self-adjoint, then for any eigenfunction y(x), the associated
eigenvalue λ must satisfy:

∫ₐᵇ y(x) · L[y(x)] dx = -λ · ∫ₐᵇ w(x) · [y(x)]² dx

Since both integrals are real, it follows that λ must be real. This property is essential for physical
systems where imaginary eigenvalues would correspond to non-physical modes.

Property 3: Completeness

Description: The set of eigenfunctions forms a complete basis in the space of square-
integrable functions over [a, b] with respect to the weight w(x).

This means that any sufficiently smooth and well-behaved function f(x) defined on [a, b]
can be expanded as an infinite series of the form:

f(x) = ∑ₙ cₙ · yₙ(x)

where the coefficients are computed via the weighted inner product:

cₙ = ∫ₐᵇ w(x) · f(x) · yₙ(x) dx

Completeness allows us to represent arbitrary functions using the Sturm–Liouville

eigenfunctions, which is particularly useful in solving partial differential equations via separation
of variables.

Property 4: Role of the Weight Function

Description: The weight function w(x) defines the inner product space and governs the
orthogonality of the eigenfunctions.

In Sturm–Liouville theory, the inner product of two functions u(x) and v(x) is defined as:

⟨u, v⟩ = ∫ₐ ᵇ w(x) · u(x) · v(x) dx

This weighted inner product ensures that the eigenfunctions are orthogonal in a way that is
intrinsic to the differential operator. The weight function arises naturally from the form of the
equation and cannot be arbitrarily chosen.
Property 5: Effect of Boundary Conditions

Description: The boundary conditions directly determine the form of eigenfunctions and
the spectrum of eigenvalues.

Consider typical boundary conditions:

● Dirichlet: y(a) = 0, y(b) = 0

● Neumann: y′(a) = 0, y′(b) = 0
● Robin (mixed): a₁ · y(a) + a₂ · y′(a) = 0

The specific type of boundary condition affects which eigenfunctions are admissible. Only
solutions that satisfy the boundary conditions are accepted, and this imposes quantization on
the spectrum of eigenvalues (i.e., discrete values λₙ ).

Property 6: Regular vs. Singular Sturm–Liouville Problem

Description: The classification into regular or singular problems depends on the behavior of
the functions p(x), q(x), w(x) and the interval [a, b].

A. Regular Sturm–Liouville Problem:

Occurs when:

● The interval [a, b] is finite,

● p(x) > 0 and continuous on [a, b],
● q(x) and w(x) > 0 are continuous on [a, b].

These ensure well-posedness with a discrete and real spectrum of eigenvalues.

B. Singular Sturm–Liouville Problem:

Occurs if:

● The interval is infinite (e.g., [0, ∞)),

● Or p(x), q(x), or w(x) are singular (e.g., have discontinuities or diverge) at
endpoints.

Singular problems often arise in physical systems with infinite domains or singularities, such as
Bessel’s or Legendre’s equations.
2.8 Applications
The Sturm–Liouville theory is a cornerstone in solving boundary value problems (BVPs)
involving linear differential equations, especially in mathematical physics and engineering. Its
applications extend to the solution of partial differential equations (PDEs) using separation
of variables, Fourier series expansions, and eigenfunction expansions.

1. Solution of Partial Differential Equations (PDEs)

One of the most powerful applications of Sturm–Liouville theory is in the separation of

variables technique used to solve PDEs such as:

● Heat Equation
● Wave Equation
● Laplace Equation

In such solutions, the spatial part of the PDE often reduces to a Sturm–Liouville problem.

Example: Heat Equation in 1D

Consider the 1D heat equation:

∂u/∂t = α² ∂²u/∂x²

Using separation of variables:

Let u(x, t) = X(x) T(t)

This leads to:

T'(t)/α²T(t) = X''(x)/X(x) = −λ
X''(x) + λX(x) = 0

with boundary conditions such as X(0) = X(L) = 0

This is a classic Sturm–Liouville problem, and the eigenfunctions Xₙ(x) form a Fourier sine
series basis to express the solution.

2. Fourier Series and Eigenfunction Expansions

The set of eigenfunctions from a regular Sturm–Liouville problem forms an orthogonal basis
for a function space with respect to a weight function w(x). Hence, any suitable function f(x)
can be expanded as:

f(x) = Σₙ cₙ yₙ(x)

where yₙ(x) are eigenfunctions and:

cₙ = ⟨f, yₙ⟩ / ⟨yₙ, yₙ⟩

= ∫ₐᵇ f(x) yₙ(x) w(x) dx / ∫ ₐᵇ yₙ²(x) w(x) dx

This eigenfunction expansion is fundamental in physics and engineering.

3. Quantum Mechanics

In quantum mechanics, the time-independent Schrödinger equation is of Sturm–Liouville

form:

−(ħ² / 2m) d²ψ/dx² + V(x)ψ = Eψ

which can be rewritten as:

d/dx (p(x) dψ/dx) + [λw(x) − q(x)] ψ = 0

The solutions ψ(x) are eigenfunctions, and E are eigenvalues representing quantized energy
levels.

Thus, Sturm–Liouville theory provides the framework for:

● Quantum particle in a box

● Harmonic oscillator
● Hydrogen atom
and more.

4. Vibrating Strings and Membranes

The displacement of a vibrating string fixed at both ends satisfies the wave equation:

∂²u/∂t² = c² ∂²u/∂x²

Applying separation of variables reduces the spatial component to:

X''(x) + λX(x) = 0
X(0) = X(L) = 0

This Sturm–Liouville problem leads to sine eigenfunctions, allowing solutions to be written as

superpositions of vibrational modes:

u(x, t) = Σ Aₙ sin(nπx/L) cos(nπct/L)

Each eigenvalue corresponds to a natural frequency of vibration.

5. Heat Conduction and Diffusion Problems

Sturm–Liouville theory is essential in solving non-steady-state conduction problems in rods

or plates. When using separation of variables, one encounters a spatial eigenvalue problem of
the Sturm–Liouville type.

These problems are often subject to non-homogeneous or inhomogeneous boundary

conditions, which determine the nature of the spectrum and eigenfunctions.

6. Signal Processing and Image Reconstruction

The orthogonal properties of Sturm–Liouville eigenfunctions are used in Fourier transforms,

wavelet analysis, and other orthogonal expansions in signal and image processing.

They provide a basis for decomposing signals into frequency components, enhancing
applications in:

● Audio compression
● Image de-noising
● Feature extraction

7. Structural Mechanics and Elasticity

In the study of beams and plates (Euler–Bernoulli beam theory), solutions to BVPs involve
eigenvalue problems of the Sturm–Liouville type.

For example, the deflection y(x) of a uniform beam under vibration satisfies:

EI y''''(x) = ρA ∂²y/∂t²
which after separation of variables leads to an ODE in Sturm–Liouville form.

2.9 MCQs
1. A Sturm–Liouville problem is considered regular if:
A) p(x) or r(x) becomes infinite in the interval
B) The interval is infinite
C) The coefficients p(x), q(x), r(x) are continuous and p(x), r(x) > 0 in a finite interval
D) The eigenfunctions are orthonormal

2. The weight function r(x) in the Sturm–Liouville problem is used:

A) To normalize the solution
B) To define orthogonality of eigenfunctions
C) To make boundary conditions homogeneous
D) To eliminate singularities

3. Eigenfunctions of a Sturm–Liouville problem are orthogonal with respect to:

A) The standard dot product
B) No known function
C) A weight function r(x)
D) The boundary values

4. Which of the following is true about the eigenvalues in a regular Sturm–Liouville problem?
A) All eigenvalues are imaginary
B) There are only a finite number of eigenvalues
C) All eigenvalues are real and form an infinite sequence
D) They depend on the choice of q(x) only

5. The completeness property of eigenfunctions implies that:

A) Only periodic functions can be represented
B) A function can be represented as a series of eigenfunctions
C) Boundary conditions must be Dirichlet
D) Orthogonality fails

6. The differential operator in a Sturm–Liouville problem is:

A) Non-linear and non-self-adjoint
B) Linear but not self-adjoint
C) Linear and self-adjoint
D) A second-order partial differential operator

7. The general Sturm–Liouville form is:

A) d²y/dx² + q(x)y = λy
B) d/dx [p(x) dy/dx] + [λr(x) − q(x)] y = 0
C) d²y/dx² = −λy
D) p(x)y'' + q(x)y = λy

8. In the inner product defined for a Sturm–Liouville problem, which expression is used?
A) ∫ₐᵇ y₁(x) y₂(x) dx
B) ∫ₐᵇ y₁′(x) y₂′(x) dx
C) ∫ₐᵇ r(x) y₁(x) y₂(x) dx
D) ∫ₐᵇ p(x) y₁″(x) y₂(x) dx

9. The effect of boundary conditions on the Sturm–Liouville problem is to:

A) Ensure orthogonality of eigenfunctions
B) Ensure uniqueness of eigenfunctions only
C) Determine the specific form of the eigenfunctions and eigenvalues
D) Affect only the weight function

10. A singular Sturm–Liouville problem occurs when:

A) p(x) or r(x) is discontinuous or infinite at the interval endpoints
B) Eigenvalues are non-real
C) The domain is symmetric
D) The weight function is zero

Answer Key

1. C
2. B
3. C
4. C
5. B
6. C
7. B
8. C
9. C
10. A

2.10 Solved Illustrations

Example 1: Verifying Sturm–Liouville Form

Problem:
Show that the differential equation
 d/dx [ (1 − x²) dy/dx ] + λy = 0
is a Sturm–Liouville problem on the interval −1 < x < 1.

Solution:

A Sturm–Liouville equation has the standard form:

d/dx [ p(x) dy/dx ] + [ λr(x) − q(x) ] y = 0

Here,
p(x) = 1 − x²
q(x) = 0
r(x) = 1

Also, the domain is finite: −1 < x < 1.

Check continuity:

● p(x) = 1 − x² is continuous on (−1, 1) and > 0

● r(x) = 1 is constant and > 0
● q(x) = 0 is continuous

Hence, this is a regular Sturm–Liouville problem.

Example 2: Orthogonality of Eigenfunctions

Problem:
Let yₙ(x) and yₘ(x) be eigenfunctions of a regular Sturm–Liouville problem:
d/dx [ p(x) dy/dx ] + [ λr(x) − q(x) ] y = 0
Show that if λₙ ≠ λₘ, then:
∫ₐᵇ r(x) yₙ(x) yₘ(x) dx = 0

Solution:

Let L[y] = −(d/dx)[p(x) dy/dx] + q(x)y

Since the operator is self-adjoint, for λₙ ≠ λₙ:

L[yₙ] = λₙ r(x) yₙ
L[yₙ] = λₙ r(x) yₙ

Multiply first by yₙ, second by yₙ, and subtract:

yₙ L[yₙ] − yₙ L[yₙ] = (λₙ − λₙ) r(x) yₙ yₙ

Integrate over [a, b]:

 ∫ₐᵇ [yₘ L[yₙ] − yₙ L[yₘ]] dx = (λₙ − λₘ) ∫ₐᵇ r(x) yₙ(x) yₘ(x) dx

The left side becomes zero due to self-adjointness and boundary conditions.

⇒(λₙ − λₘ) ∫ₐ ᵇ r(x) yₙ (x) yₘ(x) dx = 0

⇒∫ₐ ᵇ r(x) yₙ (x) yₘ(x) dx = 0

Hence, eigenfunctions corresponding to distinct eigenvalues are orthogonal with

respect to r(x).

Example 3: Solving a Simple Sturm–Liouville Problem

Problem:
Solve the Sturm–Liouville problem:
d²y/dx² + λy = 0
on the interval [0, π], with boundary conditions y(0) = 0, y(π) = 0.

Solution:

This can be rewritten as:

p(x) = 1, q(x) = 0, r(x) = 1
This is a regular Sturm–Liouville problem.

The differential equation is:

d²y/dx² + λy = 0

Case 1: λ = 0
y″ = 0 ⇒ y(x) = A + Bx
Apply BCs:
y(0) = A = 0, y(π) = A + Bπ = 0 ⇒ B = 0
⇒ y(x) ≡ 0 ⇒ Trivial Solution

Case 2: λ < 0
Let λ = −μ² ⇒ y″ − μ²y = 0
General solution: y(x) = A e^{μx} + B e^{−μx}
Apply BCs ⇒ Only trivial solution ⇒ Reject

Case 3: λ = μ² > 0
Let λ = μ² ⇒ y″ + μ²y = 0
General solution: y(x) = A sin(μx) + B cos(μx)

Apply BCs:
y(0) = B = 0
 y(π) = A sin(μπ) = 0
For non-trivial solution ⇒ sin(μπ) = 0 ⇒ μπ = nπ ⇒ μ = n

Thus,
λₙ = n², yₙ(x) = sin(n x), n = 1, 2, 3, …

Eigenvalues: λₙ = n²
Eigenfunctions: yₙ(x) = sin(n x)

Example 4: Eigenvalue Problem with Weight Function

Problem:
Solve the Sturm–Liouville problem:
d/dx [ x dy/dx ] + λx y = 0, for x ∈ (0, 1),
with boundary conditions y(0) finite, y(1) = 0.

Solution:

This is already in Sturm–Liouville form:

d/dx [ p(x) dy/dx ] + [ λr(x) − q(x) ] y = 0

Where:
p(x) = x, q(x) = 0, r(x) = x

So this is a regular Sturm–Liouville problem on (0, 1).

We rewrite the equation:

d/dx [ x dy/dx ] + λx y = 0

Let’s simplify: x y″ + y′ + λx y = 0

This is a Bessel-type equation. Try substitution:

Let y(x) = √x · u(x)
But since the equation doesn't fully reduce to known Bessel form, we use a trial solution.

Instead, observe that the eigenfunctions satisfy: x y″ + y′ + λx y = 0

This matches the Sturm–Liouville form with weight function r(x) = x.

Let us guess a solution of the form:

y(x) = xⁿ

Then,
y′ = n xⁿ⁻ ¹, y″ = n(n−1) xⁿ⁻ ²
Substitute into the equation: x [n(n−1)xⁿ⁻ ²] + n xⁿ⁻ ¹ + λx · xⁿ = 0
⇒ n(n−1)xⁿ⁻ ¹ + n xⁿ⁻ ¹ + λ xⁿ⁺ ¹ = 0
⇒ [n(n−1) + n] xⁿ⁻ ¹ + λ xⁿ⁺ ¹ = 0
⇒ n² xⁿ⁻ ¹ + λ xⁿ⁺ ¹ = 0

This cannot hold for all x unless each term separately vanishes.
So try solution by power series or numerical methods, or use known results.

Known eigenfunctions:
yₙ(x) = J₀ (αₙ x), where J₀ is the Bessel function of the first kind of order 0
Eigenvalues:
λₙ = (αₙ)², where αₙ is the nth zero of J₀ (x)

Apply BC:
y(1) = 0 ⇒ J₀ (αₙ) = 0

Thus,
Eigenfunctions: yₙ(x) = J₀ (αₙ x)
Eigenvalues: λₙ = αₙ², where αₙ is the nth zero of J₀

Example 5: Problem with Nonconstant Coefficients

Problem:
Reduce the equation
(x² y′)′ + λ x² y = 0
into standard Sturm–Liouville form.

Solution:

We are given:
(x² y′)′ + λ x² y = 0

We identify:
p(x) = x²
r(x) = x²
q(x) = 0

So the Sturm–Liouville form is:

d/dx [ p(x) dy/dx ] + [ λ r(x) − q(x) ] y = 0
This matches directly.

Hence, this is a Sturm–Liouville problem with:

● p(x) = x²
 ● q(x) = 0
● r(x) = x²

This is regular in the domain x ∈ (0, ∞) or any finite interval (a, b) with a > 0.

Let us try to solve:

(x² y′)′ + λ x² y = 0

Differentiate:
d/dx [x² y′] = 2x y′ + x² y″
So, the equation becomes:
x² y″ + 2x y′ + λ x² y = 0

Divide by x²:
y″ + (2/x) y′ + λ y = 0

This is a second-order linear ODE. Try solution:

y(x) = xⁿ

Then y′ = n xⁿ⁻ ¹, y″ = n(n−1) xⁿ⁻ ²

Substitute:
n(n−1)xⁿ⁻ ² + (2/x) n xⁿ⁻ ¹ + λ xⁿ = 0
⇒ [n(n−1) + 2n] xⁿ⁻ ² + λ xⁿ = 0
⇒ [n² + n] xⁿ⁻ ² + λ xⁿ = 0

Again, this does not satisfy the equation unless λ = 0 or x is a specific value.
So we note the standard form and try series or numerical solutions.

Conclusion:

● The equation is already in Sturm–Liouville form with

p(x) = x², q(x) = 0, r(x) = x²
● The normalized version:
y″ + (2/x) y′ + λ y = 0
This can be used for series or special function solutions.

2.11 Unsolved Problems

Problem 1:
Show that the differential operator
L[y] = d/dx(p(x) dy/dx) + λw(x)y
is self-adjoint on the interval under suitable boundary conditions, where p(x) > 0 and w(x) > 0.
Problem 2:
Find the eigenvalues and eigenfunctions of the Sturm–Liouville problem:
d/dx[(1 - x²) dy/dx] + λy = 0, for x ∈ (−1, 1),
with boundary conditions y(−1) = y(1) = 0.

Problem 3:
Determine whether the following Sturm–Liouville problem is regular or singular:
d/dx(x dy/dx) + λx y = 0, for x ∈ (0, 1),
with boundary conditions y(0) finite, y(1) = 0.

Problem 4:
Prove that eigenfunctions corresponding to distinct eigenvalues of a regular Sturm–Liouville
problem are orthogonal with respect to the weight function w(x).

Problem 5:
Convert the following boundary value problem into a Sturm–Liouville form:
y″ + λy = 0, for x ∈ (0, π),
with boundary conditions y(0) = 0, y′(π) = 0.

Problem 6:
Given the Sturm–Liouville form:
(sin x y′)′ + λ sin x y = 0, for x ∈ (0, π),
determine the appropriate orthogonality relation for the eigenfunctions.

Problem 7:
Let L[y] = −(p(x) y′)′ + q(x) y. Show that the inner product
⟨f, g⟩ = ∫ₐ ᵇ f(x)g(x)w(x) dx
ensures that L is symmetric under certain conditions on p(x), q(x), and w(x).

Problem 8:
Suppose φₙ(x) are eigenfunctions of a Sturm–Liouville problem with weight function w(x).
Prove that the set {φₙ(x)} is complete in the space of square-integrable functions over with
respect to w(x).

Problem 9:
Given a nonhomogeneous Sturm–Liouville equation:
(p(x)y′)′ + [λw(x) + q(x)]y = f(x)
with appropriate boundary conditions, describe the general strategy using eigenfunction
expansion to find the solution.

Problem 10:
For the problem:
d/dx(x² dy/dx) + λx² y = 0, for x ∈ (1, e),
with boundary conditions y(1) = 0, y(e) = 0,
determine whether the problem is regular and what is the weight function.

2.12 Summary
Sturm–Liouville theory provides a powerful and unifying approach for solving a wide range of
problems in mathematical physics. Its properties of orthogonality, completeness, and real
eigenvalues underpin many techniques for solving PDEs and modeling physical systems.

2.13 Key Terms

● Sturm–Liouville Problem
● Self-Adjoint Operator
● Eigenvalue
● Weight Function
● Orthogonality
● Completeness
● Regular/Singular Problems
● Boundary Conditions
● Eigenfunction Expansion
● Separation of Variables
2.14 References
● Boyce, W. E., & DiPrima, R. C. (2017). Elementary Differential Equations and Boundary
Value Problems. Wiley.
● Haberman, R. (2012). Applied Partial Differential Equations. Pearson.
● Zettl, A. (2005). Sturm-Liouville Theory. American Mathematical Society.
● Tenenbaum, M., & Pollard, H. (1985). Ordinary Differential Equations. Dover.
Chapter 3: Green’s Functions

Learning Objectives

By the end of this chapter, students will be able to:

● Understand the concept and mathematical definition of Green’s functions.

● Apply Green’s functions to solve second-order linear differential equations with boundary
conditions.
● Analyze the role of the Dirac delta function in the definition of Green’s function.
● Explore applications of Green’s functions across physics and engineering domains.
● Use Green’s functions to convert differential problems into integral problems.

3.1 Introduction

Many boundary value problems (BVPs) in mathematical physics, engineering, and applied
mathematics can be modeled using linear differential equations. In particular, inhomogeneous
problems of the form:

L[y(x)] = f(x)

arise frequently, where:

● L is a linear differential operator (usually second order),

● f(x) is a given source function or forcing term,
● y(x) is the unknown function to be determined.

In these situations, Green’s functions provide a structured and graceful method for addressing
these equations, particularly when homogeneous boundary conditions are applied. This
technique converts the differential problem into an integral form, thus transforming the challenge
of solving a differential equation into the task of calculating an integral. Green's function can be
viewed as the impulse response of the operator L, much like how systems react to a Dirac delta
function input in control theory or signal processing.

3.2 Mathematical Definition

Let L be a linear differential operator acting on a suitable function space over an interval [a, b],
and suppose we want to solve:

L[y(x)] = f(x), for x ∈ [a, b]

subject to homogeneous boundary conditions. Then the Green’s function, denoted by G(x, ξ),
is defined as the function that satisfies:

L[G(x, ξ)] = δ(x − ξ), for x ∈ [a, b]

with G(x, ξ) satisfying the same boundary conditions in the variable x as y(x).

Here, δ(x − ξ) is the Dirac delta function, and ξ is a fixed point in [a, b]. Conceptually, G(x, ξ)
captures the response of the system described by L to a unit impulse applied at x = ξ.

3.3 Dirac Delta Function

The Dirac delta function δ(x − ξ) is not a function in the classical sense, but a generalized
function or distribution, with the defining property:

∫ₐᵇ f(x) δ(x − ξ) dx = f(ξ), if ξ ∈ [a, b]

and δ(x − ξ) = 0 for x ≠ ξ.

In differential equations, the delta function is interpreted as a concentrated source or impulse

at x = ξ, and its action is only understood under the integral sign.

In the Green’s function equation:

L[G(x, ξ)] = δ(x − ξ)

the function G(x, ξ) must be constructed so that the effect of L acting on G produces this
idealized impulse at x = ξ.

3.4 Construction of Green’s Function for Second-Order Operators

Let us now specialize to second-order linear differential operators, which are most common
in applications. Consider:

L[y] = d²y/dx² + p(x) dy/dx + q(x) y

To find the Green’s function G(x, ξ), we must solve:

d²G/dx² + p(x) dG/dx + q(x) G = δ(x − ξ)

subject to boundary conditions in x, such as:

G(a, ξ) = 0 and G(b, ξ) = 0, or similar homogeneous conditions.

The key idea is to:

● Solve L[G] = 0 separately on (a, ξ) and (ξ, b),

● Enforce continuity of G(x, ξ) at x = ξ,
● Impose a jump condition on the first derivative dG/dx at x = ξ to accommodate the
delta function.

This leads to a piecewise-defined function G(x, ξ), which is smooth except for a discontinuity in
its derivative at x = ξ.

3.5 Theoretical Significance

Green’s functions:

● Convert differential equations into integral forms.

● Prove operator invertibility.
● Are foundational in eigenfunction expansions and Green's identities.
● Appear in fields such as quantum mechanics, heat conduction, and vibrations.

3.6 Applications of Green’s Functions

1. Solving Linear Differential Equations with Nonhomogeneous Terms

Consider a second-order linear differential equation:

L[y(x)] = f(x),
where L is a linear differential operator and f(x) is a known forcing function.

Using Green’s function G(x, ξ), the solution y(x) is given by the integral:

y(x) = ∫ G(x, ξ) f(ξ) dξ,

where the integration is carried over the domain of the problem.

Here:

● G(x, ξ) accounts for the effect at x due to a unit impulse applied at ξ.

● This method transforms a differential equation into an integral equation.

2. Electrostatics: Poisson's Equation

In electrostatics, the potential φ(x) due to a charge distribution ρ(x) satisfies the Poisson
equation:

∇ ²φ(x) = -ρ(x)/ε₀

The Green’s function G(x, ξ) satisfies:

∇ ²G(x, ξ) = -δ(x - ξ)

Then the potential is:

φ(x) = ∫ G(x, ξ) ρ(ξ)/ε₀ dξ

In 3D space, the Green’s function for free space is:

G(x, ξ) = 1 / (4π‖x - ξ‖)

Hence,

φ(x) = (1 / 4πε₀ ) ∫ ρ(ξ) / ‖x - ξ‖ dξ

This represents the superposition of point charges in continuous space.

3. Heat Conduction

The 1D heat equation is:

∂u/∂t = α ∂²u/∂x²

With initial and boundary conditions, the solution can be expressed as:

u(x, t) = ∫ G(x, ξ, t) u₀ (ξ) dξ

where:

● G(x, ξ, t) is the Green’s function for the heat equation.

● u₀ (ξ) is the initial temperature distribution.

Example for infinite domain:

G(x, ξ, t) = 1 / √(4παt) * exp[ - (x - ξ)² / (4αt) ]

This reflects how a heat pulse at ξ spreads over time.

4. Vibrating String and Wave Equation

The 1D wave equation:

∂²u/∂t² = c² ∂²u/∂x²

Solution with Green’s function approach:

u(x, t) = ∫∫ G(x, t; ξ, τ) f(ξ, τ) dξ dτ

where:

● G(x, t; ξ, τ) satisfies the homogeneous wave equation with an impulse at position ξ and
time τ.
● f(ξ, τ) represents the external force acting on the string.

5. Quantum Mechanics: Time-Independent Schrödinger Equation

In one dimension, the time-independent Schrödinger equation is:

-ħ² / (2m) * d²ψ/dx² + V(x)ψ = Eψ

This can be rearranged into:

(d²/dx² + k²)ψ(x) = (2m/ħ²)V(x)ψ(x)

Here, k² = 2mE / ħ².

The Green’s function for the operator (d²/dx² + k²) is used to express the wavefunction:

ψ(x) = ∫ G(x, ξ) V(ξ)ψ(ξ) dξ

This is the Lippmann-Schwinger equation, commonly used in quantum scattering theory.

6. Elasticity and Structural Mechanics

In structural analysis, Green’s functions can be used to compute the deflection y(x) of a beam
subjected to a distributed load q(x), modeled by:

EI d⁴y/dx⁴ = q(x)

Using the Green’s function G(x, ξ) satisfying:

 EI d⁴G/dx⁴ = δ(x - ξ),

the deflection is:

y(x) = ∫ G(x, ξ) q(ξ) dξ

This enables the superposition of effects of distributed loads on structures like beams and
plates.

7. Boundary Control and Inverse Problems

In control theory and inverse problems, Green’s functions help identify:

● How boundary inputs (e.g., heat or force applied) influence the system.
● The location or magnitude of unknown sources within a domain.

For example, if the temperature distribution u(x) is measured, and the Green’s function is
known, one may attempt to recover f(x) in:

u(x) = ∫ G(x, ξ) f(ξ) dξ

This is the essence of inverse problems, which are ill-posed but important in geophysics,
medical imaging, etc.

Check Your Progress: MCQs

1. What is the defining property of a Green’s function G(x, ξ) for a linear differential operator L?
A) L(G(x, ξ)) = 1
B) L(G(x, ξ)) = δ(x − ξ)
C) L(G(x, ξ)) = 0
D) G(x, ξ) = δ(x + ξ)

2. Which of the following best describes the role of the Dirac delta function in the Green’s
function approach?
A) It smooths out the solution.
B) It replaces the boundary condition.
C) It acts as a source term concentrated at a point.
D) It defines the eigenvalue of the system.

3. For the problem L(y) = f(x) with suitable boundary conditions, the solution using Green’s
function is given by:
A) y(x) = ∫ G(x, ξ) dξ
B) y(x) = ∫ f(ξ) dξ
C) y(x) = ∫ G(ξ, x) f(ξ) dξ
D) y(x) = f(x) × G(x, x)

4. Which condition must a Green’s function satisfy at the point x = ξ?

A) Continuity of the function and its second derivative
B) Discontinuity in the function itself
C) Continuity of the function and a jump in its first derivative
D) Jump in both the function and its derivative

5. In a self-adjoint operator L, the Green’s function satisfies:

A) G(x, ξ) = G(ξ, x)
B) G(x, ξ) = −G(ξ, x)
C) G(x, ξ) = 0
D) G(x, ξ) = 1

6. Which of the following boundary conditions allows construction of a symmetric Green’s

function?
A) Dirichlet conditions only
B) Neumann conditions only
C) Either Dirichlet or Neumann if L is self-adjoint
D) No boundary condition guarantees symmetry

7. What is the advantage of using Green’s function in solving linear inhomogeneous differential
equations?
A) It avoids the need to solve the homogeneous part.
B) It transforms the problem into a boundary value problem.
C) It allows expressing the solution as an integral involving the forcing function.
D) It removes all need for boundary conditions.

8. In the context of partial differential equations, Green’s functions are used to solve:
A) Only time-dependent equations
B) Only Laplace’s equation
C) Linear PDEs with boundary conditions
D) Nonlinear differential equations

9. The solution to a nonhomogeneous equation using Green’s function is valid if:

A) The operator is nonlinear
B) The domain is unbounded
C) The Green’s function satisfies the same boundary conditions
D) The Green’s function is arbitrary

10. If G(x, ξ) is the Green’s function for the operator L, then the solution y(x) to L(y) = f(x) is:
A) y(x) = G(x, x) × f(x)
B) y(x) = ∫ f(ξ) × G(x, ξ) dξ
C) y(x) = ∫ δ(x − ξ) × f(ξ) dξ
D) y(x) = ∫ G(x, ξ) dx

Answer Key

1. B
2. C
3. C
4. C
5. A
6. C
7. C
8. C
9. C
10. B

Key Terms

● Green’s function
● Dirac delta function
● Linear operator
● Boundary value problem
● Impulse response
● Jump condition
● Integral representation

Solved Illustrations

Illustration 1: Green’s Function for a Simple Boundary Value Problem

Problem:

Find the Green's function for the boundary value problem:

d²y/dx² = f(x),
with boundary conditions y(0) = y(1) = 0.

Solution:
This is a second-order linear differential equation with homogeneous boundary conditions.

We define the Green's function G(x, ξ) such that:

d²G(x, ξ)/dx² = δ(x − ξ),

G(0, ξ) = G(1, ξ) = 0.

We split the solution domain into two intervals:

● For x < ξ
● For x > ξ

Let: G(x, ξ) = A(ξ)⋅ x + B(ξ), for 0 ≤ x < ξ

G(x, ξ) = C(ξ)⋅ x + D(ξ), for ξ < x ≤ 1

We impose the boundary conditions: G(0, ξ) = 0 ⇒ B(ξ) = 0

G(1, ξ) = 0 ⇒ C(ξ)⋅ 1 + D(ξ) = 0 ⇒ D(ξ) = −C(ξ)

Continuity at x = ξ:
A(ξ)⋅ ξ = C(ξ)⋅ ξ + D(ξ)

Discontinuity in derivative due to δ-function:

dG/dx (x = ξ⁺ ) − dG/dx (x = ξ⁻ ) = 1

This gives: C(ξ) − A(ξ) = 1

Solving the two equations:

From continuity: A(ξ)⋅ ξ = C(ξ)⋅ ξ + D(ξ)

⇒ A(ξ)⋅ ξ = C(ξ)⋅ ξ − C(ξ)
⇒ A(ξ)⋅ ξ = C(ξ)(ξ − 1)

Now: A(ξ) = C(ξ) + 1

Substitute into the equation: (C(ξ) + 1)⋅ ξ = C(ξ)(ξ − 1)

⇒ C(ξ)⋅ ξ + ξ = C(ξ)⋅ (ξ − 1)

⇒ C(ξ)⋅ ξ + ξ = C(ξ)⋅ ξ − C(ξ)

⇒ ξ = −C(ξ) ⇒ C(ξ) = −ξ
Then, A(ξ) = −ξ + 1

Now we have:

For x < ξ: G(x, ξ) = (1 − ξ)⋅ x

For x > ξ: G(x, ξ) = −ξ⋅ x + ξ

Illustration 2: Solving a BVP Using Green’s Function

Problem:

Solve the boundary value problem:

d²y/dx² = x,
y(0) = y(1) = 0

using Green's function.

Solution:

From Illustration 1, we already found the Green’s function:

G(x, ξ) =
(1 − ξ)⋅ x for 0 ≤ x < ξ
ξ⋅ (1 − x) for ξ < x ≤ 1

Now, the solution to the BVP is:

y(x) = ∫₀ ¹ G(x, ξ)⋅ f(ξ) dξ

Here, f(ξ) = ξ

So:

y(x) = ∫₀ ˣ (1 − ξ)⋅ x⋅ ξ dξ + ∫ₓ¹ ξ⋅ (1 − x)⋅ ξ dξ

= x⋅ ∫₀ ˣ (1 − ξ)⋅ ξ dξ + (1 − x)⋅ ∫ₓ¹ ξ² dξ

Compute the integrals:

∫₀ ˣ (1 − ξ)⋅ ξ dξ
= ∫₀ ˣ (ξ − ξ²) dξ
= [ξ²/2 − ξ³/3] from 0 to x
= x²/2 − x³/3

∫ₓ¹ ξ² dξ
= [ξ³/3] from x to 1
= (1/3) − x³/3

Now substitute:

y(x) = x⋅ (x²/2 − x³/3) + (1 − x)⋅ ((1/3) − x³/3)

Simplify:

y(x) = x³/2 − x⁴/3 + (1 − x)⋅ (1/3 − x³/3)

= x³/2 − x⁴/3 + (1/3 − x³/3 − x/3 + x⁴/3)

= x³/2 − x⁴/3 + 1/3 − x³/3 − x/3 + x⁴/3

Now combine like terms:

−x⁴/3 + x⁴/3 = 0
x³/2 − x³/3 = x³⋅ (1/2 − 1/3) = x³/6

Final result:

y(x) = x³/6 − x/3 + 1/3

Solved Illustration 3

Problem: Solve the boundary value problem using Green’s function:

d²y/dx² = sin(x), with y(0) = 0, y(π) = 0.

Solution:

We need to find y(x) such that: L(y) = d²y/dx² = sin(x), with y(0) = y(π) = 0.

The Green’s function G(x, ξ) for the operator L = d²/dx² with Dirichlet BCs (y = 0 at both ends)
over [0, π] is given by:

G(x, ξ) =
(ξ(π − x))/π for x < ξ
(x(π − ξ))/π for x ≥ ξ

Then, the solution is: y(x) = ∫₀ ^π G(x, ξ) sin(ξ) dξ

Let’s split the integral at x:

y(x) = ∫₀ ^x (ξ(π − x)/π) sin(ξ) dξ + ∫ₓ^π (x(π − ξ)/π) sin(ξ) dξ

This can be evaluated using integration by parts or numerically for a particular value of x. For
example, at x = π/2:
 First integral: ∫₀ ^{π/2} (ξ(π − π/2)/π) sin(ξ) dξ = (π/2π) ∫₀ ^{π/2} ξ sin(ξ) dξ
= (1/2) ∫₀ ^{π/2} ξ sin(ξ) dξ
= (1/2) [−ξ cos(ξ) + ∫ cos(ξ) dξ]₀ ^{π/2}
= (1/2)[−(π/2) cos(π/2) + sin(π/2) − (0 − 0)] = (1/2)(0 + 1) = 0.5

Second integral: Similar procedure gives ≈ 0.5

So, y(π/2) ≈ 0.5 + 0.5 = 1.0

Solved Illustration 4

Problem: Solve using Green’s function:

d²y/dx² = 1, with y(0) = y(1) = 0

Solution:

We want the solution to: d²y/dx² = 1, on [0, 1], with y(0) = y(1) = 0

Green’s function for this BVP is: G(x, ξ) =

ξ(1 − x) for x < ξ
x(1 − ξ) for x ≥ ξ

Then the solution is: y(x) = ∫₀ ^1 G(x, ξ) ⋅ 1 dξ

= ∫₀ ^x ξ(1 − x) dξ + ∫ₓ^1 x(1 − ξ) dξ

First integral:
∫₀ ^x ξ(1 − x) dξ = (1 − x) ∫₀ ^x ξ dξ = (1 − x)(x²/2)

Second integral:
∫ₓ^1 x(1 − ξ) dξ = x ∫ₓ^1 (1 − ξ) dξ = x[(1 − ξ)²/2]ₓ^1 = x[(1 − x)²/2]

So, y(x) = (1 − x)(x²/2) + x(1 − x)²/2

= (1 − x)x²/2 + x(1 − x)²/2

Factor: y(x) = (1/2)[x²(1 − x) + x(1 − x)²]

= (1/2)x(1 − x)[x + (1 − x)]
= (1/2)x(1 − x)(1)

Final Answer: y(x) = (1/2) x (1 − x)

Self-Assessment Questions
Instructions: Solve the following boundary value problems using Green’s function methods.
Express the general form of the Green’s function and the integral solution. Wherever applicable,
perform explicit integration.

Problem 1

Solve using Green’s function:

d²y/dx² = x, with y(0) = 0, y(1) = 0

Problem 2

Use Green’s function to solve:

d²y/dx² = eˣ, with y(0) = y(π) = 0

Problem 3

Solve the BVP using Green’s function:

d²y/dx² = cos(x), with y(−π) = 0, y(π) = 0

Problem 4

Find the Green’s function and solve:

d²y/dx² = 2x + 3, with y(0) = 0, y(2) = 0

Problem 5

Solve the BVP using Green’s function:

d²y/dx² = δ(x − 1), with y(0) = y(2) = 0
(Here, δ is the Dirac delta function)

Problem 6
Consider the non-homogeneous equation:
d²y/dx² − y = x, with y(0) = y(1) = 0
Construct the corresponding Green’s function and evaluate the solution.

Problem 7

Let L(y) = d²y/dx² on [0, 1] with y(0) = 0, y′(1) = 0.

Find the Green’s function corresponding to L and use it to solve L(y) = sin(πx).

Problem 8

Find the solution using Green’s function for:

d²y/dx² = ln(x), with y(1) = y(2) = 0
(Note: ln(x) is defined and integrable over [1, 2])

Problem 9

Let L(y) = d²y/dx² with Neumann boundary conditions:

y′(0) = 0, y′(1) = 0.
Construct the Green’s function and solve L(y) = x²

Problem 10

Solve using Green’s function:

d²y/dx² + y = cos(x), with y(0) = y(π/2) = 0
Determine the appropriate Green’s function and the integral solution.

References (APA Style)

● Green, G. (1828). An Essay on the Application of Mathematical Analysis to the Theories

of Electricity and Magnetism.
● Arfken, G. B., & Weber, H. J. (2005). Mathematical Methods for Physicists (6th ed.).
Elsevier.
● Haberman, R. (2013). Applied Partial Differential Equations with Fourier Series and
Boundary Value Problems (5th ed.). Pearson.
● Strauss, W. A. (2007). Partial Differential Equations: An Introduction (2nd ed.). Wiley.
Chapter 4: Applications of Boundary Value Problems

Learning Objectives

By the end of this chapter, learners will be able to:

● Distinguish BVPs from IVPs in real-world modeling.

● Understand the mathematical form of BVPs and why they model equilibrium systems.
● Identify and classify boundary conditions (Dirichlet, Neumann, Mixed).
● Recognize BVP applications in heat transfer, vibration analysis, electrostatics, fluid flow,
and quantum mechanics.
● Interpret eigenvalue problems arising in BVP contexts.

4.1 What Makes BVPs Special?

In contrast to initial value problems (IVPs), which set conditions at a single location (typically the
beginning of the range), boundary value problems (BVPs) define conditions at two or more
separate points. This approach more effectively represents real-world systems in areas such as
physics, engineering, and applied mathematics.

Example:

To model a string fixed at both ends, you do not need how the string started vibrating (an initial
condition), but you must ensure that its displacement is zero at both fixed ends — this is a
boundary condition.

4.2 Mathematical Form of BVP

A standard second-order linear BVP is written as:

(1)
d/dx ( p(x) * dy/dx ) + q(x) * y = λ * w(x) * y

with boundary conditions such as:

(2)
a₁ * y(a) + a₂ * y′(a) = 0
b₁ * y(b) + b₂ * y′(b) = 0

Where:

● p(x), q(x), w(x) are given coefficient functions.

 ● λ is a spectral parameter (often an eigenvalue).
● y(x) is the unknown function to be solved.
● a and b define the endpoints of the domain [a, b].

4.3 Why Are BVPs Widely Applicable?

BVPs naturally arise in equilibrium or steady-state problems in physics and engineering, where
the conditions they create are of more interest than their temporal evolution.

Examples of common uses include steady-state thermal distribution

● Elastic beams bending statically

● Potential of electricity in a conductor
● Patterns of vibration in mechanical systems
● Particle quantum states contained within a box

4.4 Key Features of BVPs

(a) Spatial Domain Control

Boundary Value Problems (BVPs) are phenomena such as temperature along a rod or voltage
along a wire that are located over a spatial range and do not change over time.

b) Physics Realism

Boundary conditions, such as fixed endpoints, insulated boundaries, or surfaces of zero flux, are
often used to represent physical constraints.

c) Eigenvalue-Related Issues

Many BVPs become eigenvalue issues, especially in wave mechanics or vibration systems
where certain frequencies (eigenvalues) produce nontrivial solutions.

4.5 General Application Example

Steady-State Heat Equation in One Dimension

Consider a rod of length L with thermal conductivity k and heat sources/sinks along its length.
The governing BVP is:
d²T/dx² + f(x)/k = 0, 0<x<L

with boundary conditions like:

T(0) = T₀ , T(L) = T_L

Here, T(x) is the temperature, and f(x) is a heat source function.

This represents a BVP because the temperature is fixed at both ends, and we seek the
equilibrium profile T(x).

4.6 Types of Boundary Conditions

a) Dirichlet Boundary Conditions

Specify the value of the function at the boundary.

ₙ T(0) = T₀ , T(L) = T_L

b) Neumann Boundary Conditions

Specify the derivative (i.e., the flux) at the boundary.

ₙ dT/dx(0) = q₀ , dT/dx(L) = q_L

c) Mixed Boundary Conditions

Combine both value and derivative.

ₙ T(0) + α * dT/dx(0) = C

These reflect different physical setups — fixed temperature, insulated boundaries, or convective
surfaces.

4.7 Importance in Eigenfunction Expansions

In many cases, the solutions to BVPs form an orthogonal set of eigenfunctions. These
eigenfunctions are used to expand arbitrary input functions (e.g., heat distribution, force
distribution) into a Fourier-like series.

This is foundational in:

● Solving partial differential equations (PDEs)

● Performing modal analysis in engineering
● Understanding quantum mechanical wave functions
Applications of BVPs
4.8 Heat Conduction in a Rod

Equation: d²u/dx² = 0
BCs: u(0) = T₁ , u(L) = T₂
Solution: u(x) = [(T₂ - T₁ )/L]·x + T₁

4.9 Vibrating String

Wave equation: ∂²u/∂t² = c²·∂²u/∂x²

BCs: u(0, t) = u(L, t) = 0

Using separation of variables, spatial ODE becomes: d²X/dx² + λ·X = 0, X(0) = X(L) = 0

Eigenvalues: λₙ = (n·π/L)²
Eigenfunctions: Xₙ(x) = sin(n·π·x / L)

Full solution:
u(x, t) = Σ [Aₙ·cos(n·π·c·t/L) + Bₙ·sin(n·π·c·t/L)] · sin(n·π·x / L)

4.10 Electrostatics

Laplace Equation: ∂²u/∂x² + ∂²u/∂y² = 0

BCs on edges → solution via separation of variables:

● u(x, y) = Σ Aₙ·sin(n·π·x/a)·sinh(n·π·y/a)

4.11 Quantum Mechanics: Particle in a Box

Equation: d²ψ/dx² + k²·ψ = 0, ψ(0) = ψ(L) = 0

Eigenvalues: Eₙ = (n²·π²·ħ²)/(2m·L²)
Eigenfunctions: ψₙ(x) = √(2/L)·sin(n·π·x/L)
4.12 Fluid Flow: Poiseuille Flow

Equation: d²u/dy² = -G/μ

BCs: u(0) = u(h) = 0

Solution: u(y) = (G/2μ)·y·(h − y)

Summary

● BVPs are essential in modeling steady-state and spatially constrained phenomena.

● Solutions are often framed via eigenvalue problems and orthogonal expansions.
● Physical interpretations span thermodynamics, vibrations, quantum states, and flow
mechanics.

Key Terms

● Dirichlet Condition
● Neumann Condition
● Mixed Boundary Condition
● Eigenvalue
● Poiseuille Flow
● Sturm–Liouville System

Self-Assessment Questions

1. In the context of a vibrating string, the solution to the wave equation subject to fixed boundary
conditions is typically expressed using:
A) Fourier series
B) Laplace transform
C) Green’s function
D) Runge-Kutta method

2. Which type of boundary value problem is most associated with steady-state heat conduction?
A) Hyperbolic BVP
B) Elliptic BVP
C) Parabolic BVP
D) Dynamic BVP
3. In the application of BVPs to electrostatics, the potential field u(x, y) often satisfies:
A) ∂u/∂t = α ∂²u/∂x²
B) ∇ ²u = 0
C) d²u/dt² + ω²u = 0
D) ∇ ²u = −ρ

4. The method of separation of variables is commonly used in solving which kind of BVPs?
A) Linear algebraic equations
B) Ordinary differential equations only
C) Partial differential equations with homogeneous boundary conditions
D) Nonlinear PDEs with random forcing

5. The use of Fourier series in boundary value problems requires:

A) Initial conditions only
B) Boundary conditions to be discontinuous
C) Periodicity or domain boundedness
D) Use of complex numbers

6. Which of the following is not a typical application of boundary value problems?

A) Quantum mechanics
B) Structural engineering
C) Random walk simulation
D) Fluid dynamics

7. The classical heat equation under Dirichlet boundary conditions leads to which form of
general solution?
A) Polynomial solutions
B) Gaussian function
C) Trigonometric series in space and exponential in time
D) Bessel function solutions

8. When using eigenfunction expansions to solve BVPs, the completeness property ensures:
A) Uniqueness of the solution
B) Finite number of eigenvalues
C) That any function satisfying boundary conditions can be represented
D) Orthogonality of the differential operator

9. Which mathematical tool helps solve inhomogeneous BVPs using impulse responses?
A) Newton-Raphson method
B) Laplace transform
C) Green’s function
D) Euler’s method

10. In fluid dynamics, the Navier–Stokes equation becomes a BVP when:

A) Time derivative is removed for steady state
B) Pressure is ignored
C) Flow is turbulent
D) Flow is uniform in all directions

Answer Key

1. A
2. B
3. B
4. C
5. C
6. C
7. C
8. C
9. C
10. A

Solved Illustrations

Solved Illustration 1: Heat Equation in a Rod with Dirichlet Boundary

Conditions

Problem:
Solve the one-dimensional heat equation
∂u/∂t = α ∂²u/∂x², for 0 < x < L, t > 0
with boundary conditions:
u(0, t) = 0, u(L, t) = 0
and initial condition:
u(x, 0) = sin(πx/L)

Solution:
We solve using separation of variables.
Assume u(x, t) = X(x) T(t)

Substitute into the PDE:

(1/T) dT/dt = α (1/X) d²X/dx² = −λ

So we have two ODEs:

d²X/dx² + λX = 0
dT/dt + αλT = 0
Apply boundary conditions on X(x):
X(0) = 0, X(L) = 0
⇒ The eigenvalues are: λₙ = (nπ/L)², n = 1, 2, 3, ...

Eigenfunctions: Xₙ(x) = sin(nπx/L)

Time solution: Tₙ(t) = e^(−α(nπ/L)² t)

Now apply initial condition:

u(x, 0) = sin(πx/L)
⇒ Only the first term is present in the Fourier expansion

Final solution:
u(x, t) = sin(πx/L) × e^(−α(π/L)² t)

Solved Illustration 2: Vibration of a Stretched String

Problem:
Solve the wave equation for a string of length L fixed at both ends:
∂²u/∂t² = c² ∂²u/∂x², 0 < x < L, t > 0
with boundary conditions:
u(0, t) = 0, u(L, t) = 0
and initial conditions:
u(x, 0) = sin(2πx/L), ∂u/∂t(x, 0) = 0

Solution:
Use separation of variables:
Let u(x, t) = X(x) T(t)

Substitute:
(1/T) d²T/dt² = c² (1/X) d²X/dx² = −λ

Get two ODEs:

d²X/dx² + λX = 0
d²T/dt² + c²λ T = 0

Apply boundary conditions:

X(0) = 0, X(L) = 0
⇒ λₙ = (nπ/L)², Xₙ(x) = sin(nπx/L)

So for n = 2 (due to initial condition),

u(x, t) = sin(2πx/L) × cos(2πct/L)
 Solved Illustration 3: Steady-State Temperature Distribution in a Plate

Problem:
Solve Laplace’s equation ∇ ²u = 0 on a square plate 0 < x < 1, 0 < y < 1, with:

● u(0, y) = 0, u(1, y) = 0
● u(x, 0) = 0
● u(x, 1) = 100

Solution:
Assume u(x, y) = X(x) Y(y)

Substitute into Laplace’s equation:

d²X/dx² / X = −d²Y/dy² / Y = −λ

Two ODEs:
X″ + λX = 0 with X(0) = 0, X(1) = 0
⇒ Xₙ(x) = sin(nπx), λ = (nπ)²

Y″ − λY = 0
⇒ Yₙ(y) = Aₙ sinh(nπy)

General solution:
u(x, y) = Σ Aₙ sin(nπx) sinh(nπy)

Apply boundary condition u(x, 1) = 100:

Use Fourier sine series to expand 100 in terms of sin(nπx)

This leads to:

Aₙ = (4 / (nπ)) for odd n
⇒ u(x, y) = Σ_{n odd} (4 / (nπ)) sin(nπx) sinh(nπy) / sinh(nπ)

Would you like a few more illustrations or shall I proceed to the unsolved problem set?

Certainly. Here are three more solved illustrations for the chapter Applications of Boundary
Value Problems, with all equations in Unicode plain text format for direct copy-pasting into
Google Docs:

Solved Illustration 4: Deflection of a Beam (Euler–Bernoulli Equation)

Problem:
Find the deflection y(x) of a beam of length L simply supported at both ends, subjected to a
uniform load q₀ .
The governing equation is:
d⁴y/dx⁴ = q₀ / EI
with boundary conditions:
y(0) = 0, y(L) = 0
y″(0) = 0, y″(L) = 0
Here, E = modulus of elasticity, I = moment of inertia.

Solution:
Integrate the differential equation:
d⁴y/dx⁴ = q₀ / EI

First integration:
d³y/dx³ = (q₀ /EI)x + C₁
Second:
d²y/dx² = (q₀ /2EI)x² + C₁ x + C₂
Third:
dy/dx = (q₀ /6EI)x³ + (C₁ /2)x² + C₂ x + C₃
Fourth:
y(x) = (q₀ /24EI)x⁴ + (C₁ /6)x³ + (C₂ /2)x² + C₃ x + C₄

Apply boundary conditions:

● y(0) = 0 ⇒ C₄ = 0
● y(L) = 0 ⇒ (q₀ /24EI)L⁴ + (C₁ /6)L³ + (C₂ /2)L² + C₃ L = 0
● y″(0) = 0 ⇒ C₂ = 0
● y″(L) = 0 ⇒ (q₀ /2EI)L² + C₁ L = 0 ⇒ C₁ = −(q₀ /2EI)L

Now plug back values to find y(x):

Final solution:
y(x) = (q₀ / 24EI) [x⁴ − 2Lx³ + L²x²]

Solved Illustration 5: Electric Potential in a Rectangular Domain

Problem:
Solve Laplace’s equation ∇ ²φ = 0 in a rectangle 0 < x < a, 0 < y < b, with:

● φ(0, y) = φ(a, y) = φ(x, 0) = 0

● φ(x, b) = V₀

Solution:
Let φ(x, y) = X(x) Y(y)
Then:
X″/X = −Y″/Y = −λ
So,
X″ + λX = 0 with X(0) = 0, X(a) = 0
⇒ λₙ = (nπ/a)², Xₙ(x) = sin(nπx/a)

Y″ − λY = 0
⇒ Yₙ(y) = Aₙ sinh(nπy/a)

General solution:
φ(x, y) = Σ Aₙ sin(nπx/a) sinh(nπy/a)

Apply boundary condition φ(x, b) = V₀:

Expand V₀ as a sine series:
Aₙ = (2V₀ / (nπ sinh(nπb/a))) × ∫₀ᵃ sin(nπx/a) dx = (4V₀ / (nπ sinh(nπb/a))) for odd n

Final solution:
φ(x, y) = Σ_{n odd} (4V₀ / (nπ sinh(nπb/a))) × sin(nπx/a) sinh(nπy/a)

Solved Illustration 6: Eigenvalue Problem in Sturm–Liouville Form

Problem:
Solve the Sturm–Liouville problem:
(d/dx)[p(x) dy/dx] + [λ w(x)] y = 0, 0<x<π
with: p(x) = 1, w(x) = 1, and boundary conditions y(0) = y(π) = 0.

Solution:
Substitute p(x) = w(x) = 1 into the equation:
d²y/dx² + λy = 0

This is a standard eigenvalue problem.

Apply boundary conditions:
y(0) = 0, y(π) = 0

The general solution of the ODE:

y(x) = A sin(√λ x) + B cos(√λ x)

Apply y(0) = 0 ⇒ B = 0
Apply y(π) = 0 ⇒ sin(√λ π) = 0 ⇒ √λ π = nπ ⇒ λ = n², n = 1, 2, 3, ...

Eigenvalues: λₙ = n²
Eigenfunctions: yₙ(x) = sin(n x)
Unsolved Problem Set: Applications of Boundary Value
Problems

Problem 1: Heat Conduction in a Rod

A rod of length L is perfectly insulated along its lateral surface. The temperature u(x, t) in the rod
satisfies the heat equation:
∂u/∂t = α² ∂²u/∂x², 0 < x < L, t > 0

Initial condition: u(x, 0) = sin(πx / L)

Boundary conditions: u(0, t) = u(L, t) = 0

Find the temperature distribution u(x, t) for t > 0.

Problem 2: Beam Deflection with Fixed-Free Ends

A cantilever beam of length L is fixed at x = 0 and free at x = L. The beam is subject to a point
load P at the free end. The governing equation is:
EI d⁴y/dx⁴ = 0

Boundary conditions:
y(0) = 0, dy/dx(0) = 0
d²y/dx²(L) = 0, d³y/dx³(L) = −P/EI

Find the equation of the deflection curve y(x).

Problem 3: Electrostatic Potential in a Square Region

Solve Laplace’s equation in a square region 0 < x < 1, 0 < y < 1:

∇ ²φ = 0

Boundary conditions:
φ(x, 0) = 0, φ(x, 1) = 0
φ(0, y) = 0, φ(1, y) = V₀

Find the steady-state potential φ(x, y).

Problem 4: Eigenvalue Problem with Variable Coefficients

Solve the eigenvalue problem:

d/dx [ (1 + x²) dy/dx ] + λ y = 0, 0<x<1
with boundary conditions:
y(0) = 0, y(1) = 0

Find the eigenvalues λ and corresponding eigenfunctions y(x).

Problem 5: Green’s Function for a Differential Operator

Let L[y] = −d²y/dx² on the interval 0 < x < 1, with boundary conditions y(0) = y(1) = 0.
Find the Green’s function G(x, ξ) that satisfies:
L[G(x, ξ)] = δ(x − ξ)

Problem 6: Population Dynamics Model

Consider a population density u(x, t) on a region 0 < x < L governed by:

∂u/∂t = D ∂²u/∂x² + ru(1 − u/K)

Assume Neumann boundary conditions (no flux):

∂u/∂x(0, t) = ∂u/∂x(L, t) = 0

Formulate this as a boundary value problem and discuss its qualitative behavior.

Problem 7: Fluid Flow in a Channel

The velocity profile u(y) of fluid flowing between two plates is governed by:
μ d²u/dy² = dp/dx

where dp/dx is a constant pressure gradient and μ is the viscosity.

Boundary conditions: u(0) = u(h) = 0

Find the velocity distribution u(y).

Problem 8: Vibrating Membrane – Circular Domain

Consider the radial part R(r) of the vibrating circular membrane equation:
r² d²R/dr² + r dR/dr + (λr² − n²) R = 0

with R(a) = 0, where a is the radius of the membrane.

Identify the type of differential equation and the nature of its solutions.

References (APA Style)

● Boyce, W. E., & DiPrima, R. C. (2017). Elementary Differential Equations and Boundary
Value Problems (11th ed.). Wiley.
● Haberman, R. (2012). Applied Partial Differential Equations with Fourier Series and
Boundary Value Problems (5th ed.). Pearson.
● Strauss, W. A. (2007). Partial Differential Equations: An Introduction (2nd ed.). Wiley.

End of Chapter 4