Chapter 3

Series solutions and special functions

3.1 Power series solutions

So far we have been understanding solutions of differential equations whose coefficients and in-
homogeneous terms are all continuous in an open interval and describing methods to find their
solution. But in many situations, the coefficients of a differential equation and the inhomogeneous
terms are all in fact analytic in an open interval. In such situations, the solutions to ordinary differ-
ential equations are also analytic. Hence we may assume that they have a power series expansion.
Thus solving the differential equation boils down to figuring out the coefficients in the power series
expansion. This method is called the power series method and we shall understand that next. The
method is quite simple, and it entails the following.

Theorem 3.1.1. Suppose we have the following differential equation,

h(x)y →→ + p(x)y → + q(x)y = r(x) , (3.1)

where h(x), p(x), q(x) and r(x) are all analytic in a neighborhood of a point x0 and such that
h(x0 ) ̸= 0. Then every solution of the differential equation (3.1) is analytic in the an open disc
containing x0 .

Since we shall be working with the assumption that h(x0 ) ̸= 0, we shall suppose that h is never
zero in the domain of interest and divide it out. This leads to an equation of the form

y →→ + p(x)y → + q(x) = r(x) . (3.2)

3.1.1 Method of solution

The method of solution for solving an equation of the form (3.2) is quite straightforward. From
Theorem 3.1.1, if it is given that h(x)p(x), q(x) and r(x) are all analytic around a point x0 with h
non-vanishing, then we may suppose that a solution y is of the form
∞
!
y(x) = an (x − x0 )n .
n=0

This gives us
∞
!
y → (x) = nan (x − x0 )n−1
n=1

41
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

and
∞
!
y →→ (x) = n(n − 1)an (x − x0 )n−2
n=2
Substituting this back into the differential equation and equating for powers of x allows us to solve
for ai ’s. This in turn will allow us to obtain y.

3.2 Legendre polynomials

The Legendre differential equation associated to a parameter n is given by

(1 − x2 )y →→ − 2xy → + n(n + 1)y = 0 . (3.3)

Using the power series method, we guess that the solution is of the form
∞
!
y(x) = an xn .
n=0

We may then derive the recursive relation,

(n − s)(n + s + 1)
as+2 = − as
(s + 2)(s + 1)
for s ! 0.
can show that the solution is of the form

y(x) = a0 y1 (x) + a1 y2 (x),

where
n(n + 1) 2 (n − 2)n(n + 1)(n + 3) 4
y1 (x) = 1 − x + x − ...,
2! 4!
and
(n − 1)(n + 2) 3 (n − 3)(n − 1)(n + 2)(n + 4) 5
y2 (x) = x − x + x − ....
3! 5!
Generally, the series defining y1 , y2 are convergent for |x| < 1, but for some parameters n, we
have some additional features. When the parameter of the Legendre equation n is a non-negative
integer, the series defining either y1 or y2 reduces to a polynomial and they are called Legendre
polynomials.
The general form of the n-th Legendre polynomial is
n
$2%
! (2n − 2m)!
Pn (x) = (−1)m xn−2m . (3.4)
2n m!(n − m)!(n − 2m)!
m=0

3.3 Bessel functions

3.3.1 Frobenius method
Definition 3.3.1. A regular point of the differential equation (3.1) is a point where the functions
h(x), p(x), q(x), r(x) are all analytic and h(x0 ) ̸= 0. Any point that is not a regular point is called
a singluar point.

42
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

When a differential equation has a singular point in the domain, we cannot use the power series
method to solve the problem. We have to modify it. This modification is known as the Frobenius
method. To demonstrate this method, we try to solve the following differential equation.

b(x) → c(x)
y →→ + y + 2 y=0.
x x
Note that the point x = 0 is a singular point of the equation. We guess a solution of the form
∞
!
y(x) = xr an xn
n=0

where r ∈ C and a0 ̸= 0. Since we have assumed that the coefficients are analytic in a neighborhood
around zero, they are expressed as a series. Suppose their series expansions are given by
∞
!
b(x) = b i x i = b 0 + b 1 x + b2 x 2 + . . .
i=0

and
∞
!
c(x) = c j xj = c 0 + c 1 x + c 2 x2 + . . . .
j=0

We calculate the derivatives as

∞
!
y → (x) = an (n + r)xn+r−1
n=0

and
∞
!
→→
y (x) = an (n + r)(n + r − 1)xn+r−2 .
n=0
Substituting this back into
b(x) → c(x)
y →→ + y + 2 y = 0,
x x
we get the first term as
∞
!
y →→ = an (n + r)(n + r − 1)xn+r−2 ,
n=0
the second term
" ∞
#" ∞
# ∞ !
∞
b(x) → 1 ! ! !
i n+r−1
y = bi x an (n + r)x = bi an (n + r)xi+n+r−2
x x
i=0 n=0 i=0 n=0

and the third term as

⎛ ⎞" #
∞ ∞ ∞ !
∞
c(x) 1 ⎝! j⎠
!
n+r
!
y = c j x a n x = cj an xj+n+r−2 .
x2 x2
j=0 n=0 j=0 n=0

Write all terms as powers of xn+r−2 . For the second term let k = n + i:
∞
" k #
! !
bi ak−i (k − i + r) xk+r−2
k=0 i=0

43
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

For the third term let k = n + j:

⎛ ⎞
∞
! k
!
⎝ cj ak−j ⎠ xk+r−2
k=0 j=0

Thus the equation becomes

" k # ⎛ ⎞
∞
! ∞
! ! ∞
! k
!
ak (k + r)(k + r − 1)xk+r−2 + bi ak−i (k − i + r) xk+r−2 + ⎝ cj ak−j ⎠ xk+r−2 = 0
k=0 k=0 i=0 k=0 j=0

The lowest power of x occurs when k = 0:

a0 r(r − 1)xr−2 + b0 a0 rxr−2 + c0 a0 xr−2 .

Thus the coefficient of xr−2 is

a0 [r(r − 1) + b0 r + c0 ] .
Since a0 ̸= 0, we obtain the indicial equation

r(r − 1) + b0 r + c0 = 0. (3.5)

The roots of this equation determine the possible values of r for the Frobenius solution. Let the
roots of the indicial equation be r1 and r2 , with r1 ≥ r2 .

Case 1: Distinct Roots Not Differing by an Integer

If r1 − r2 is not an integer, then there exist two linearly independent Frobenius solutions
∞
!
y1 (x) = an xn+r1 ,
n=0

and
∞
!
y2 (x) = An xn+r2 .
n=0
Thus the general solution is
y(x) = C1 y1 (x) + C2 y2 (x).

Case 2: Repeated Root

If the indicial equation has a repeated root r1 = r2 = r, then one solution is
∞
!
y1 (x) = an xn+r .
n=0

A second linearly independent solution has the form

∞
!
y2 (x) = y1 (x) ln x + bn xn+r .
n=0

Thus
y(x) = C1 y1 (x) + C2 y2 (x).

44
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Case 3: Roots Differ by a Positive Integer

If r1 − r2 = N where N is a positive integer, then one solution is
∞
!
y1 (x) = an xn+r1 .
n=0

The second solution may either be

∞
!
y2 (x) = An xn+r2 ,
n=0

or may involve a logarithmic term of the form

∞
!
y2 (x) = C y1 (x) ln x + bn xn+r2 .
n=0

Thus the general solution is

y(x) = C1 y1 (x) + C2 y2 (x).

3.3.2 The Bessel equation

An important second–order differential equation that arises in many physical applications is the
Bessel equation. It appears when solving partial differential equations such as the Laplace
equation, heat equation, and wave equation in cylindrical or polar coordinates.
The Bessel equation of order ω is given by

x2 y →→ + xy → + (x2 − ω 2 )y = 0, (3.6)

where ω ! 0 is a constant parameter called the order of the Bessel equation. The point x = 0 is a
singular point of this differential equation. Therefore, solutions near x = 0 can be obtained using
the Frobenius method. If we divide the equation by x2 , we obtain
( )
1 → ω2
→→
y + y + 1 − 2 y = 0, (3.7)
x x
which shows explicitly that the coefficients have the form required for the application of the Frobe-
nius method. Applying the Frobenius method leads to the indicial equation

r(r − 1) + r − ω 2 = 0,

which simplifies to
r2 − ω 2 = 0. (3.8)

3.3.3 Bessel functions and properties

Thus the roots of the indicial equation are

r1 = ω, r2 = −ω.

The corresponding series solutions define the Bessel functions of the first kind

Jν (x) and J−ν (x).

45
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Consider the Bessel equation of order ω ! 0. We shall obtain a series expansion for Jν (x) using
the Frobenius method. We start with the guess solution,
∞
!
y(x) = ak xk+r , a0 ̸= 0.
k=0

Computing derivatives, we have

∞
!
→
y (x) = ak (k + r)xk+r−1 ,
k=0

and
∞
!
y →→ (x) = ak (k + r)(k + r − 1)xk+r−2 .
k=0

Substituting this back, we get

x2 y →→ + xy → + (x2 − ω 2 )y = 0.
First term:
∞
!
2 →→
x y = ak (k + r)(k + r − 1)xk+r .
k=0

Second term:
∞
!
xy → = ak (k + r)xk+r .
k=0

Third term:
∞
! ∞
!
2 2 k+r+2 2
(x − ω )y = ak x −ω ak xk+r .
k=0 k=0

Combining the first three sums gives

∞
! ∞
!
* 2
+ k+r
ak (k + r)(k + r − 1) + (k + r) − ω x + ak xk+r+2 = 0.
k=0 k=0

Simplifying,
∞
! ∞
!
* 2 2
+ k+r
ak (k + r) − ω x + ak xk+r+2 = 0.
k=0 k=0

Shift the index in the second sum by letting k → k − 2:

∞
! ∞
!
* 2 2
+ k+r
ak (k + r) − ω x + ak−2 xk+r = 0.
k=0 k=2

From the lowest power (k = 0):

a0 (r2 − ω 2 ) = 0.
This is the indicial equation as before. For k ≥ 2, we get the recursive relation
* +
ak (k + r)2 − ω 2 + ak−2 = 0. (3.9)

46
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Thus
ak−2
ak = − .
(k + r)2 − ω 2
For simplicity, let us suppose that
Choose r = n. Then

(k + n)2 − n2 = k(k + 2n).

Hence the recurrence relation becomes
ak−2
ak = − .
k(k + 2n)
Only even powers appear in the series. Let k = 2m:

(−1)m a0
a2m = .
22m m!(n + m)!
Thus the solution becomes

(−1)m , x -2m
∞
!
y(x) = a0 xn .
m!(n + m)! 2
m=0

Choosing the normalization

1
a0 = ,
2n n!
we obtain the Bessel function of the first kind

(−1)m , x -2m+n
∞
!
Jn (x) = . (3.10)
m!(n + m)! 2
m=0

This function is called the Bessel function of the first kind of order n.

3.3.4 The Gamma Function

In the series expansion of Bessel functions of non–integral order, expressions involving factorials
appear in which the argument is not necessarily an integer. Since the factorial n! is defined only
for non–negative integers, it is necessary to introduce a function that extends the factorial to more
general values. This extension is provided by the Gamma function. The Gamma function is
defined for α > 0 by the improper integral
. ∞
!(α) = xα−1 e−x dx.
0

This function satisfies the fundamental recurrence relation

!(α + 1) = α !(α).

If α = n is a positive integer, repeated application of this relation gives

!(n + 1) = n!.

47
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Thus the Gamma function extends the factorial function to non–integer arguments. For this reason
it is often regarded as the natural generalization of the factorial.
Using the Gamma function, factorial expressions can be written in the form

n! = !(n + 1),

and more generally expressions such as (ω + m)! that appear in the series solutions of differential
equations can be interpreted as
(ω + m)! = !(ω + m + 1),
even when ω is not an integer.
This notation allows the series expansion of the Bessel function of order ω to be written in a
unified form for both integer and non–integer values of ω.

3.3.5 Bessel Function of the First Kind for Non–Integral Order

Using the Gamma function, the series solution of the Bessel equation can be written in a form that
is valid for non–integral values of the order ω. Recall that the Bessel equation of order ω is

x2 y →→ + xy → + (x2 − ω 2 )y = 0.

Applying the Frobenius method at the regular singular point x = 0, we assume a solution of the
form
!∞
y(x) = ak xk+r , a0 ̸= 0.
k=0

Substitution into the differential equation gives the indicial equation

r2 − ω 2 = 0,

so that
r1 = ω, r2 = −ω.
Consider the solution corresponding to r = ω. The recurrence relation for the coefficients obtained
from the Frobenius method is
ak−2
ak = − .
(k + ω)2 − ω 2
Since
(k + ω)2 − ω 2 = k(k + 2ω),
we obtain
ak−2
ak = − .
k(k + 2ω)
As before, only even powers occur in the series. Writing k = 2m and solving the recurrence relation
gives
(−1)m a0
a2m = 2m .
2 m!!(m + ω + 1)
Hence the corresponding series solution becomes
∞
! (−1)m , x -2m
y(x) = a0 xν .
m!!(m + ω + 1) 2
m=0

48
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Choosing the normalization

1
a0 = ,
2ν !(ω + 1)
we obtain the Bessel function of the first kind of order ω
∞
! (−1)m , x -2m+ν
Jν (x) = .
m!!(m + ω + 1) 2
m=0

This series is valid for any real or complex value of ω. When ω = n is a non–negative integer, the
Gamma function satisfies
!(n + m + 1) = (n + m)!,
and the series reduces to the previously obtained expansion for the Bessel function of integer order.

3.3.6 Properties of Bessel functions

The Bessel functions of the first kind satisfy several useful recurrence relations and derivative
identities.

Basic Recurrence Relations

The Bessel functions satisfy the relations
2ω
Jν−1 (x) + Jν+1 (x) = Jν (x),
x
and
Jν−1 (x) − Jν+1 (x) = 2Jν→ (x).
These formulas relate Bessel functions of neighboring orders.

Derivative Identities
From the recurrence relations we obtain the derivative formulas
1
Jν→ (x) = (Jν−1 (x) − Jν+1 (x)) ,
2
and
d
(xν Jν (x)) = xν Jν−1 (x),
dx
and
d / −ν 0
x Jν (x) = −x−ν Jν+1 (x).
dx

Alternative Recurrence Forms

The recurrence relations can also be written in the form
2ω
Jν+1 (x) = Jν (x) − Jν−1 (x),
x
and
2ω
Jν−1 (x) = Jν (x) − Jν+1 (x).
x

49
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Derivative in Terms of Neighboring Orders

Another useful derivative formula is
d ω
Jν (x) = Jν (x) − Jν+1 (x),
dx x
or equivalently
d ω
Jν (x) = − Jν (x) + Jν−1 (x).
dx x
We also have the following derivative relations, An important identity is
d
(xν Jν (x)) = xν Jν−1 (x).
dx
More generally,
d / −ν 0
x Jν (x) = −x−ν Jν+1 (x).
dx
These relations are frequently used when solving boundary value problems and when manipulating
expressions involving Bessel functions.

3.3.7 General solution of the Bessel equation and the Bessel functions of the
second kind
If ω is not an integer, then two Frobenius solutions are obtained,
∞
! (−1)m , x -2m+ν
Jν (x) = ,
m!!(m + ω + 1) 2
m=0

and
∞
! (−1)m , x -2m−ν
J−ν (x) = .
m!!(m − ω + 1) 2
m=0

When ω is not an integer, these two functions are linearly independent and form a fundamental
set of solutions of the Bessel equation (3.6). However, suppose ω = n, where n is an integer. Using
properties of the Gamma function, it can be shown that

J−n (x) = (−1)n Jn (x).

Thus J−n (x) and Jn (x) are linearly dependent.

Therefore, when ω = n is an integer, the Frobenius method yields only one linearly independent
solution. Since a second–order differential equation requires two linearly independent solutions to
form the general solution, we must seek another solution that is independent of Jn (x).
This leads to the introduction of the Bessel functions of the second kind, denoted by Yn (x),
which together with Jn (x) form a fundamental set of solutions of the Bessel equation for integer
order.
y(x) = C1 Jn (x) + C2 Yn (x).

50
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Solved Exercises
Exercise 3.1. Solve the differential equation y → = y.

1 Since the
Solution. coefficients are all analytic, we may employ the power series method. Suppose
n is the power series expansion of y around the point x = 0. Then from the given
y(x) = ∞ n=0 a n x 0
condition, we have
!∞ !∞
nan xn−1 = am xm .
n=1 m=1
In particular, a1 = a0 , a2 = a1
2 = a0
2 etc. More generally , we can say that an = n! .
a0
Thus
∞
! a0
y(x) = xn = a0 ex .
n!
n=0

Exercise 3.2. Solve the differential equation y →→ + y = 0 using the power series method.
Solution. Once again, we may suppose that y is a power series around the point x0 = 0 with
coefficients an . Substituting this back into the differential equation, we have
∞
! ∞
!
an n(n − 1)xn−2 = − am xm .
n=2 m=0

Hence 2a2 = −a0 , 6a3 = −a1 , 12a4 = −a2 = and so on giving us that
a0
2
( ) ( )
x2 x4 x3 x5
y(x) = a0 1 − + − . . . + a1 x − + − . . . = a0 sin(x) + a1 cos(x).
2! 4! 3! 5!
Exercise 3.3. Solve the differential equation y →→ − y = x.
Solution. We use the power series method and expand y as a power series around the point x0 = 0.
Expanding this out, we get
∞
! ∞
!
n−2
an n(n − 1)x − am xm = x.
n=2 m=0

In other words,
∞
!
(am+2 (m + 2)(m + 1) − am ) xm = x.
m=0
Equating the coefficients, we get 2a2 − a0 = 0, 6a3 − a1 = 1, 30a4 − a2 = 0, 20a5 − a3 = 0 and so
on. More generally we have the formula,
am−2 am−4 a0
am = = = ... = ,
m(m − 1) m(m − 1)(m − 2)(m − 3) m!
if m is even. If m ! 3 is odd, then
a1 + 1
am = .
m!
Plugging this back we get
∞
! ∞
!
x2m x2m+1
y(x) = a0 + (a1 + 1) + a1 x.
(2m)! (2m + 1)!
m=0 m=1

51
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Rewriting the above we get that

( ) ( )
ex + e−x ex − e−x
y(x) = a0 + (a1 + 1) −x .
2 2
√
Exercise 3.4. Show that !(1/2) = π.
Solution. From the definition of the Gamma function,
.∞
/ 0
! 21 = x−1/2 e−x dx.
0

Make the substitution x = t2 , so that dx = 2t dt. Then

. ∞
/10 2
! 2 = (t2 )−1/2 e−t (2t dt).
0

Since (t2 )−1/2 = 1/t, this simplifies to

. ∞ .∞
/ 0 −t2 2
! 12 = 2 e dt = e−t dt.
0
−∞

Let . ∞
2
I= e−t dt.
−∞
Then (. ∞ ) (. ∞ ) .
2 2 2 +y 2 )
I2 = e−x dx e−y dy = e−(x dx dy.
−∞ −∞ R2
Switching to polar coordinates (r, θ),
. 2π . ∞
2
I2 = e−r r dr dθ.
0 0

Evaluating the integrals, . ∞

2
2
I = 2π re−r dr.
0
Let u = r2 . Then du = 2r dr, giving
I 2 = π.
Hence √
I= π.
Since . √
∞
−t2 π
e dt = ,
0 2
we obtain / 0 √
! 12 = π.
Exercise 3.5. Show that 2
2
J1/2 (x) = sin(x)
πx

52
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Solution. The Bessel function of the first kind has the series expansion
∞
! (−1)m , x -2m+ν
Jν (x) = .
m! !(m + ω + 1) 2
m=0
Let ω = 2.
1
Then
∞
! (−1)m , x -2m+ 1
2
J1/2 (x) = / 3
0 .
m! ! m + 2
2
m=0
Using the identity
/ 0 (2m + 1)!! √
! m + 32 = π,
2m+1
and simplifying the coefficients, the series becomes
2 ∞
2 ! (−1)m x2m+1
J1/2 (x) = .
πx (2m + 1)!
m=0
But the series
∞
! (−1)m x2m+1
(2m + 1)!
m=0
is the Taylor series for sin x. Hence
2
2
J1/2 (x) = sin x.
πx
Exercise 3.6. Suppose ω is not an integer. What is the Wronskian of Jν (x) and J−ν (x)?
Solution. Consider the second-order linear ODE
( )
→→ 1 → ω2
y + y + 1 − 2 y = 0,
x x
and let y1 , y2 be two solutions. The Wronskian is defined as
W = y1 y2→ − y2 y1→ .
Differentiate W with respect to x:
dW / 0
= y1→ y2→ + y1 y2→→ − y2→ y1→ + y2 y1→→ = y1 y2→→ − y2 y1→→ .
dx
Using the differential equation,
( )
→→ 1 → ω2
yi = − yi − 1 − 2 yi (i = 1, 2),
x x
we substitute:
dW
= y1 (−P y2→ − Qy2 ) − y2 (−P y1→ − Qy1 ).
dx
Simplifying, ( )
dW 1 → → ω2
= − (y1 y2 − y2 y1 ) − 1 − 2 (y1 y2 − y2 y1 ).
dx x x
Since y1 y2 − y2 y1 = 0, this reduces to
dW 1
= − W (x).
dx x
So W (x) = C
x, for some constant C.

53
Complex analysis and diff. eq.
CHAPTER 3. SERIES SOLUTIONS AND SPECIAL FUNCTIONS

Exercise 3.7. Suppose that ω is not an integer. Show that the zeros of the Jν Bessel functions are
simple.

Solution. Suppose for the sake of contradiction, that there is a point x0 such that J(x0 ) = J → (x0 ) =
0. But then we would have W (Jν , J−ν )(x0 ) = 0, which is not possible by the above exercise.

Exercise 3.8. Prove the interlacing properties of the zeros of the Bessel function.

Solution. Let Jν (x) be the Bessel function of the first kind, and let its positive zeros be

0 < jν,1 < jν,2 < · · · .

Since Jν (x) is continuous and differentiable, Rolle’s theorem implies that between any two
successive zeros jν,k and jν,k+1 , there exists at least one point ck such that

Jν→ (ck ) = 0.

Now use the identity

1/ 0
Jν→ (x) = Jν−1 (x) − Jν+1 (x) .
2
Thus, at x = ck ,
Jν−1 (ck ) = Jν+1 (ck ).
Next, consider the relation
2ω
Jν−1 (x) + Jν+1 (x) = Jν (x).
x
At a zero x = jν,k , we have
Jν−1 (jν,k ) = − Jν+1 (jν,k ).
Hence, as x moves from jν,k to jν,k+1 , the functions Jν−1 (x) and Jν+1 (x) change sign relative
to each other. Hence at least one of Jν−1 or Jν+1 has a zero in between jν,k and jν,k+1 . Using the
recurrence relation
2ω
Jν−1 (x) + Jν+1 (x) = Jν (x),
x
we get at these points

Jν−1 (jν,k ) = − Jν+1 (jν,k ), Jν−1 (jν,k+1 ) = − Jν+1 (jν,k+1 ).

Hence both Jν−1 and Jν+1 has to change sign in the interval (jν,k , jν,k+1 ). Since both Jν−1 (x) and
Jν+1 (x) are continuous functions, a sign change in each implies that each must vanish at least once
in the interval
(jν,k , jν,k+1 ).
Thus, both Jν−1 (x) and Jν+1 (x) have at least one zero between consecutive zeros of Jν (x). More-
over, there is exactly one zero of Jν+1 (x) and one of Jν−1 (x) lies between consecutive zeros of
Jν (x).
Thus, the zeros interlace:

jν,k < jν+1,k < jν,k+1 , jν,k < jν−1,k+1 < jν,k+1 .

54
Chapter 4

Laplace transform

4.1 Definition and properties

The Laplace transform is a powerful integral transform used to solve differential equations, espe-
cially initial value problems arising in engineering and physics. It converts differential equations
into algebraic equations, which are often easier to solve.
From a conceptual point of view, the Laplace transform converts a function of time into a
function of a new variable s, which may be interpreted as a frequency or growth–decay parameter.
In this sense, the transform changes the perspective from the time domain, where a process evolves
step by step, to a domain where the same process is described in terms of exponential modes.
Exponential functions est play a fundamental role because differentiation acts very simply on them.
The Laplace transform therefore expresses a complicated time-dependent function as a weighted
combination of these exponential behaviors.
This viewpoint explains why the Laplace transform is so useful for solving differential equations.
In the time domain, differentiation is an operation that increases complexity, but under the Laplace
transform it becomes multiplication by the variable s. As a result, a differential equation involving
derivatives of a function is converted into an algebraic equation involving the transform of that
function. Algebraic equations are generally much easier to solve, and once the solution is obtained
in the s-domain, the inverse Laplace transform converts it back into the time domain. In this
way, the Laplace transform provides a powerful method for reducing problems about change and
dynamics into problems of algebra.

Definition 4.1.1. Let f (t) be a function defined for t ≥ 0. The Laplace transform of f (t) is
defined by
. ∞
L{f (t)} = F (s) = e−st f (t) dt, (4.1)
0

provided the integral converges. Here s is a real or complex parameter, and F (s) is called the
Laplace transform of f (t).

The Laplace transform exists if the function f (t) satisfies certain conditions. In particular, if
f (t) is piecewise continuous on every finite interval [0, T ] and if there exist constants M and a such
that
|f (t)| ≤ M eat , t ≥ 0,
then the Laplace transform of f (t) exists for all s > a. A function satisfying this condition is said to
be of exponential order. The Laplace transform converts a function of the time variable t into a

57
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

function of the complex variable s. This transformation is particularly useful because differentiation
in the time domain corresponds to algebraic operations in the s-domain.

As we shall see the Laplace transform can be used to convert the problem of solving
differential equations to that of solving algebraic equations.

The Laplace transform is commonly written as F (s) = L{f (t)}. The function f (t) is called the
original function, while F (s) is called the transform function.

4.1.1 Properties of the Laplace transform

The Laplace transform satisfies several important properties that make it a powerful tool for solving
differential equations. Let
L{f (t)} = F (s), L{g(t)} = G(s).

Linearity
The Laplace transform is a linear operator. For constants a and b,

L{af (t) + bg(t)} = aL{f (t)} + bL{g(t)}.

Thus
L{af (t) + bg(t)} = aF (s) + bG(s).

First Shifting Theorem (Exponential Shift)

If
L{f (t)} = F (s),
then
L{eat f (t)} = F (s − a).

Differentiation in the Time Domain

If L{f (t)} = F (s), then
L{f → (t)} = sF (s) − f (0).

Proof. From the definition of the Laplace transform,

. ∞
→
L{f (t)} = e−st f → (t) dt.
0

We evaluate this integral using integration by parts. Let

u = e−st , dv = f → (t) dt.

Then
du = −se−st dt, v = f (t).
Applying integration by parts,
. ∞ . ∞
* +∞
e−st f → (t) dt = e−st f (t) 0 + s e−st f (t) dt.
0 0

58
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

Since f (t) is of exponential order, the term e−st f (t) tends to 0 as t → ∞. Hence
* −st +∞
e f (t) 0 = 0 − f (0) = −f (0).
Therefore . ∞
→
L{f (t)} = −f (0) + s e−st f (t) dt.
0
But . ∞
e−st f (t) dt = F (s).
0
Thus
L{f → (t)} = sF (s) − f (0).

Similarly,
L{f →→ (t)} = s2 F (s) − sf (0) − f → (0).
In general,
L{f (n) (t)} = sn F (s) − sn−1 f (0) − sn−2 f → (0) − · · · − f (n−1) (0).

Integration in the Time Domain

If L{f (t)} = F (s), then
3. t 4
F (s)
L f (τ ) dτ = . (4.2)
0 s

Proof. Let . t
g(t) = f (τ ) dτ.
0
By definition of the Laplace transform,
. ∞ . ∞ (. t )
−st −st
L{g(t)} = e g(t) dt = e f (τ ) dτ dt.
0 0 0
We now interchange the order of integration. The region of integration is 0 ≤ τ ≤ t < ∞, so
reversing the order gives
. ∞. t . ∞. ∞
−st
e f (τ ) dτ dt = e−st f (τ ) dt dτ.
0 0 0 τ
Since f (τ ) does not depend on t, it can be taken outside the inner integral:
. ∞ (. ∞ )
−st
= f (τ ) e dt dτ.
0 τ
Evaluating the inner integral,
. ∞ 5 6∞
1 e−sτ
e−st dt = − e−st = .
τ s τ s
Substituting this result, we obtain
. ∞ .
e−sτ 1 ∞ −sτ F (s)
L{g(t)} = f (τ ) dτ = e f (τ ) dτ = .
0 s s 0 s
This completes the proof.

59
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

Differentiation in the Transform Domain

If L{f (t)} = F (s), then
d
L{tf (t)} = − F (s) .
ds
More generally,
dn
L{tn f (t)} = (−1)n F (s) . (4.3)
dsn

Second Shifting Theorem (Time Shift)

The unit step function (also called the Heaviside function) is a function that models a sudden
change occurring at a specific time. It is typically denoted by u(t − a) and is defined by
7
0, t < a,
u(t − a) = (4.4)
1, t ≥ a.

Thus the function is zero before the time t = a and jumps to one at t = a. It is useful for describing
signals or forces that “turn on’’ at a particular instant.
The unit step function plays an important role in Laplace transform methods because it allows
piecewise-defined functions to be written in a compact algebraic form. If a function begins acting
only after some time a, it can often be written as u(t − a)f (t − a). If L{f (t)} = F (s), then

L{f (t − a)u(t − a)} = e−as F (s), (4.5)

where u(t − a) is the unit step function. The factor e−as represents a delay in time.

Proof. Starting from the definition of the Laplace transform,

. ∞
L{u(t − a)f (t − a)} = e−st u(t − a)f (t − a) dt.
0

Since the unit step function satisfies u(t − a) = 0 for t < a and u(t − a) = 1 for t ≥ a, the integral
becomes . ∞
e−st f (t − a) dt.
a
Now make the substitution
τ = t − a, t = τ + a, dt = dτ.
When t = a, we have τ = 0, and as t → ∞, τ → ∞. Thus
. ∞ . ∞
−st
e f (t − a) dt = e−s(τ +a) f (τ ) dτ.
a 0
Factor out the constant exponential term:
. ∞
−as
=e e−sτ f (τ ) dτ.
0

But the remaining integral is exactly the Laplace transform of f (t), namely F (s). Therefore,
L{u(t − a)f (t − a)} = e−as F (s).
This completes the proof.

60
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

Periodic functions
A function f (t) is said to be periodic with period T > 0 if

f (t + T ) = f (t)

for all t ≥ 0. Many signals and physical processes exhibit periodic behavior, so it is useful to know
how to compute their Laplace transforms efficiently.
If f (t) is periodic with period T , then its Laplace transform can be expressed in terms of a
single period of the function. The result is
. T
1
L{f (t)} = e−st f (t) dt.
1 − e−sT 0

Thus, instead of integrating over the entire interval [0, ∞), it suffices to integrate over one period
and multiply by the factor 1−e1−sT .
This formula follows by writing the Laplace transform integral as a sum over successive intervals
of length T :
. ∞ !∞ . (n+1)T
−st
e f (t) dt = e−st f (t) dt.
0 n=0 nT

Using the periodicity f (t + T ) = f (t) and performing a suitable change of variables, each integral
can be written in terms of the same integral over [0, T ]. The resulting series forms a geometric
series whose sum produces the factor 1−e1−sT .
This theorem is especially useful for functions defined piecewise over one period, such as square
waves, sawtooth waves, and other periodic signals that arise in engineering.

Convolution Theorem
The convolution of two functions is defined by
. t
(f ∗ g)(t) = f (τ )g(t − τ ) dτ. (4.6)
0

Then
L{(f ∗ g)(t)} = F (s)G(s). (4.7)

Thus convolution in the time domain corresponds to multiplication in the transform domain.

Proof. Consider the Laplace transform of the convolution:

. ∞ (. t )
−st
L{(f ∗ g)(t)} = e f (τ )g(t − τ ) dτ dt.
0 0

We interchange the order of integration. The region of integration is 0 ≤ τ ≤ t < ∞, so reversing

the order gives
. ∞. t . ∞. ∞
e−st f (τ )g(t − τ ) dτ dt = e−st f (τ )g(t − τ ) dt dτ.
0 0 0 τ

61
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

Since f (τ ) does not depend on t, it may be taken outside the inner integral:
. ∞ (. ∞ )
−st
= f (τ ) e g(t − τ ) dt dτ.
0 τ

Now make the substitution u = t − τ, t = u + τ, dt = du. When t = τ , u = 0, and when t → ∞,

u → ∞. Hence
. ∞ . ∞ . ∞
−st −s(u+τ ) −sτ
e g(t − τ ) dt = e g(u) du = e e−su g(u) du.
τ 0 0

Substituting this into the previous expression gives

. ∞ . ∞
−sτ
L{(f ∗ g)(t)} = f (τ )e dτ e−su g(u) du.
0 0

But these integrals are precisely the Laplace transforms of f (t) and g(t). This completes the
proof.

4.1.2 Table of Laplace transforms

Function f (t) Laplace Transform L{f (t)}

1
1
s
1
t
s2
n!
tn
sn+1
!(a + 1)
ta , a > 0
sa+1
1
eat
s−a
a
sin(ωt)
s2 + ω 2
s
cos(ωt)
s + ω2
2

a
sinh(at)
s − a2
2

s
cosh(at)
s − a2
2

b
eat sin(bt)
(s − a)2 + b2
s−a
eat cos(bt)
(s − a)2 + b2
e−as
u(t − a)
s

62
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

4.2 The Inverse Laplace Transform

The Laplace transform converts a function f (t) into a function F (s) by the integral
. ∞
F (s) = L{f (t)} = e−st f (t) dt.
0

In many applications we begin with the transform F (s) and wish to recover the original function
f (t). The operation that reconstructs f (t) from F (s) is called the inverse Laplace transform.
If
F (s) = L{f (t)},

then we write
f (t) = L−1 {F (s)}.

Thus the inverse Laplace transform reverses the effect of the Laplace transform and restores the
original function.

4.2.1 The Inversion Theorem

Under suitable conditions, the original function f (t) can be recovered from its Laplace transform
F (s) by the following integral formula.

Theorem 4.2.1 (Inversion theorem). Let f (t) be piecewise continuous on [0, ∞) and of exponential
order. If
F (s) = L{f (t)},

then for t > 0

γ+i∞
.
1
f (t) = est F (s) ds,
2πi
γ−i∞

where γ is a real constant chosen so that the vertical line Re(s) = γ lies to the right of all singularities
of F (s).

Although the inversion theorem provides a theoretical formula for recovering f (t), it is rarely
used directly in elementary applications. In practice, inverse Laplace transforms are usually ob-
tained by

• using tables of known transforms,

• applying partial fraction decomposition,

• using properties of the Laplace transform.

These techniques allow the inverse transform to be computed efficiently without evaluating the
Bromwich integral explicitly.

63
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

4.3 Applications to ODEs

4.3.1 Initial value problems
An initial value problem (IVP) for an ordinary differential equation consists of a differential equation
together with specified values of the unknown function and possibly some of its derivatives at
a particular point. These additional conditions are called initial conditions. While solving a
differential equation typically produces a general solution containing arbitrary constants, the role
of the initial conditions is to determine the specific values of those constants. Thus, solving an
IVP amounts to finding the precise constants in the general solution so that the resulting function
satisfies the given conditions.
For example, consider the constant–coefficient differential equation

y →→ − 3y → + 2y = 0.

The characteristic equation is (r2 − 3r + 2 = 0), which factors as ((r − 1)(r − 2) = 0). Hence the
general solution is
y(t) = C1 et + C2 e2t ,
where (C1 ) and (C2 ) are arbitrary constants. At this stage the solution represents a family of
functions, each corresponding to different choices of the constants.
Suppose we impose the initial conditions (y(0) = 1) and (y → (0) = 0). Substituting into the
general solution gives (C1 + C2 = 1), while differentiating yields (y → (t) = C1 et + 2C2 e2t ), and the
condition (y → (0) = 0) gives (C1 + 2C2 = 0). Solving this system determines the constants uniquely,
producing a single function that satisfies both the differential equation and the initial conditions.
This illustrates how an IVP selects one particular solution from the family given by the general
solution.
In general, if the differential equation is of order n, then its general solution contains n arbitrary
constants. Consequently, to determine a unique solution we must specify n initial conditions,
typically of the form
y(t0 ), y → (t0 ), y →→ (t0 ), . . . , y (n−1) (t0 ).
These conditions fix the values of the n constants appearing in the general solution. Without
enough initial conditions, the solution remains a family of functions rather than a single uniquely
determined function.

4.3.2 The subsidiary equation

The Laplace transform provides a systematic method for solving linear ordinary differential equa-
tions by converting them into algebraic equations in a new variable. When the Laplace transform
is applied to both sides of a differential equation, derivatives of the unknown function are trans-
formed into algebraic expressions involving the Laplace variable s and the initial conditions. As
a result, the original differential equation is converted into a simpler subsidiary equation in terms
of the Laplace transform of the unknown function. This algebraic equation can then be solved
using standard methods, and the desired solution of the original differential equation is obtained
by taking the inverse Laplace transform.

4.3.3 Solving ODEs

To solve an ordinary differential equation using the Laplace transform, we proceed through a
sequence of systematic steps. First, take the Laplace transform of both sides of the differential

64
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

equation. Use the standard formulas for transforms of derivatives, such as

L{y → (t)} = sY (s) − y(0), L{y →→ (t)} = s2 Y (s) − sy(0) − y → (0),

where Y (s) = L{y(t)}. Substituting these expressions converts the differential equation into an
algebraic equation in terms of Y (s) and the given initial conditions.
Next, simplify the resulting algebraic expression to obtain the subsidiary equation, which is an
equation involving Y (s) and the variable s. Solve this equation for Y (s) using ordinary algebraic
methods. If necessary, decompose the expression into simpler fractions (for example using partial
fractions) so that each term corresponds to a known Laplace transform pair.
Finally, take the inverse Laplace transform of Y (s) to obtain the solution y(t) of the origi-
nal differential equation. This step uses standard tables or known formulas for inverse Laplace
transforms. The resulting function y(t) satisfies both the differential equation and the given initial
conditions.
Example 4.3.1. Consider the non–homogeneous differential equation

y →→ + y = sin t,

with initial conditions y(0) = 0 and y → (0) = 0. We therefore apply the Laplace transform. Taking
transforms of both sides and using
1
L{y → (t)} = sY (s) − y(0), L{y →→ (t)} = s2 Y (s) − sy(0) − y → (0), L{sin t} = ,
s2 +1
we obtain
1
(s2 Y (s)) + Y (s) = .
s2 +1
Hence the subsidiary equation is
1
(s2 + 1)Y (s) = .
s2 + 1
Solving for Y (s), we get
1
Y (s) = .
(s2 + 1)2
To compute the inverse Laplace transform, we use the known transform pair
2s 1
L{t sin t} = , L{sin t − t cos t} = .
(s2 + 1)2 (s2 + 1)2

Hence
1/ 0
sin t − t cos t .
y(t) =
2
Thus the solution of the given non–homogeneous differential equation is
1
y(t) = (sin t − t cos t) ,
2
which satisfies both the differential equation and the initial conditions.

65
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

Solved Exercises
Exercise 4.1. Solve the differential equation with constant coefficients

y →→ − 3y → + 2y = 0,

together with the initial conditions y(0) = 1 and y → (0) = 0.

Solution. We apply the Laplace transform to both sides of the equation. Using the formulas

L{y → (t)} = sY (s) − y(0), L{y →→ (t)} = s2 Y (s) − sy(0) − y → (0),

we obtain
(s2 Y (s) − s − 0) − 3(sY (s) − 1) + 2Y (s) = 0.
Simplifying this expression gives

s2 Y (s) − s − 3sY (s) + 3 + 2Y (s) = 0,

and hence
(s2 − 3s + 2)Y (s) = s − 3.
This algebraic relation is the subsidiary equation. Solving for Y (s), we obtain
s−3
Y (s) = .
(s − 1)(s − 2)

We now decompose Y (s) into partial fractions:

s−3 2 1
= − .
(s − 1)(s − 2) s−1 s−2

Taking the inverse Laplace transform and using the formula

3 4
−1 1
L = eat ,
s−a

we obtain
y(t) = 2et − e2t .
This function satisfies both the differential equation and the given initial conditions.

Exercise 4.2. Evaluate the Laplace transform of the sine function.

Solution. Recall that . ∞

L{f (t)} = e−st f (t) dt.
0
Thus, . ∞
L{sin t} = e−st sin t dt.
0
To evaluate the integral, we use integration by parts. Let
. ∞
I= e−st sin t dt.
0

66
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

Choose
u = sin t, dv = e−st dt,
so that
1
du = cos t dt, v = − e−st .
s
Applying integration by parts,
5 6∞ . ∞
1 1
I = − e−st sin t + e−st cos t dt.
s 0 s 0

The boundary term vanishes, since e−st → 0 as t → ∞ for s > 0, and sin 0 = 0. Hence
.
1 ∞ −st
I= e cos t dt.
s 0

Let . ∞
J= e−st cos t dt.
0

Applying integration by parts again with

u = cos t, dv = e−st dt,

we obtain 5 6∞ . ∞
1 1
J = − e−st cos t − e−st sin t dt.
s 0 s 0

Evaluating the boundary term gives

1 1
J= − I.
s s
Substituting this expression for J into the earlier equation I = 1s J, we obtain
( )
1 1 1 1 1
I= − I = 2
− 2 I.
s s s s s

Rearranging, ( )
1 1
I 1+ 2 = .
s s2
Hence
1
I= .
s2 + 1
Therefore,
1
L{sin t} = , s > 0.
s2 +1

Exercise 4.3. Solve the non–homogeneous differential equation

y →→ − 3y → + 2y = et ,

together with the initial conditions y(0) = 0 and y → (0) = 0.

67
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

Solution. Taking the Laplace transform of both sides and using

L{y → (t)} = sY (s) − y(0), L{y →→ (t)} = s2 Y (s) − sy(0) − y → (0),

we obtain
1
(s2 Y (s)) − 3(sY (s)) + 2Y (s) = .
s−1
Hence,
1
(s2 − 3s + 2)Y (s) = .
s−1
Solving for Y (s) gives

1 1
Y (s) = = .
(s − 1)(s2 − 3s + 2) (s − 1)2 (s − 2)

We decompose this into partial fractions:

1 1 1 1
2
= − − .
(s − 1) (s − 2) s − 2 s − 1 (s − 1)2

Taking inverse Laplace transforms and using the formulas

3 4 3 4
1 1
L−1 = eat , L−1 = teat ,
s−a (s − a)2

we obtain the solution

y(t) = e2t − et − tet .
This function satisfies the given non–homogeneous differential equation and the initial conditions.

Exercise 4.4. Sovle the differential equation

y →→ + y = u(t − π),

with the initial conditions

y(0) = 0, y → (0) = 0,
where u(t − π) is the unit step (Heaviside) function.

Solution. Because the forcing term suddenly turns on at t = π, the equation is naturally described
using piecewise functions, making it inconvenient to solve by elementary methods such as unde-
termined coefficients. The Laplace transform method, however, handles such discontinuous inputs
systematically. Taking the Laplace transform of both sides and using

e−as
L{y →→ } = s2 Y (s) − sy(0) − y → (0), L{y} = Y (s), L{u(t − a)} = ,
s
we obtain
e−πs
(s2 + 1)Y (s) = .
s
Hence
e−πs
Y (s) = .
s(s2 + 1)

68
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

We first compute the inverse transform of

1
.
s(s2 + 1)

Using partial fractions,

1 1 s
= − 2 .
s(s2 + 1) s s +1
Thus 3 4
−1 1
L 2
= 1 − cos t.
s(s + 1)
Using the Heaviside shifting theorem,

L−1 {e−πs F (s)} = u(t − π)f (t − π),

we obtain
y(t) = u(t − π) [1 − cos(t − π)] .
Thus the solution is 7
0, 0 ≤ t < π,
y(t) =
1 − cos(t − π), t ≥ π.
This example illustrates how the Laplace transform naturally handles differential equations with
delayed or discontinuous forcing terms.

Exercise 4.5. Solve the differential equation

√
y →→ + y = t 2
,

with initial conditions

y(0) = 0, y → (0) = 0.
√
Solution. A forcing term such as t 2 is not a polynomial, exponential, or trigonometric function,
so it cannot be handled easily by classical methods such as undetermined coefficients. However,
the Laplace transform still applies because it works for a much wider class of functions.
Taking the Laplace transform of both sides and using

L{y →→ } = s2 Y (s) − sy(0) − y → (0), L{y} = Y (s),

together with the general formula

!(α + 1)
L{tα } = , α > −1,
sα+1
we obtain √
2 !( 2 + 1)
(s + 1)Y (s) = √ .
s 2+1

Hence √
!( 2 + 1)
Y (s) = √ .
s 2+1 (s2 + 1)

69
Complex analysis and diff. eq. CHAPTER 4. LAPLACE TRANSFORM

Thus the Laplace transform reduces the differential equation to an algebraic expression for Y (s).
The solution y(t) can then be written as the inverse Laplace transform
7 8
√ 1
y(t) = !( 2 + 1) L−1 √ .
s 2+1 (s2 + 1)

Although the inverse transform may involve special functions or convolution integrals, the Laplace
transform method still provides a systematic way to express the solution, even for forcing terms
that fall outside the scope of elementary solution techniques.

70
Chapter 5

Introduction to PDEs

In many physical, biological, and engineering systems, the state of a quantity depends on more
than one independent variable. For example, temperature in a plate depends on both spatial
coordinates and possibly time. While ordinary differential equations (ODEs) describe systems
with a single independent variable, partial differential equations (PDEs) are required when a
dependent variable depends on multiple independent variables.
In standard notation for PDEs with two independent variables, the dependent variable is written
as z = z(x, y), where x and y are the independent variables. The first–order partial derivatives are
commonly denoted by
∂z ∂z
p= , q= ,
∂x ∂y
and the second–order derivatives by
∂2z ∂2z ∂2z
r= , s= , t= .
∂x2 ∂x∂y ∂y 2
These notations are widely used when expressing PDEs compactly and when classifying equations
by their order and structure.
In the setting of two independent variables, a partial differential equation (PDE) is an
equation involving an unknown function z = z(x, y) of the variables x and y, together with some
of its partial derivatives. A PDE can often be written in the general form
F (x, y, z, p, q) = 0,
where F is a given function. The order of the PDE is the order of the highest derivative appearing
in the equation.
A solution of a PDE is a function z(x, y) that satisfies the equation when its partial derivatives
are substituted into it. In many applications, the PDE alone is not enough to determine a unique
solution; additional conditions must be specified. One common type is an initial value problem,
in which the value of the unknown function (and possibly some of its derivatives) is prescribed
along a given curve in the (x, y)–plane. These initial conditions, together with the PDE, determine
the evolution of the solution away from the curve and are essential for selecting the physically or
mathematically relevant solution among many possible ones.

5.1 Formation of PDEs

One common way of obtaining a partial differential equation (PDE) is by eliminating arbitrary
constants or arbitrary functions from a given family of relations. Suppose a relation involves a

75
Complex analysis and diff. eq. CHAPTER 5. INTRODUCTION TO PDES

dependent variable z = z(x, y) and contains one or more arbitrary constants. By taking partial
derivatives of the relation with respect to the independent variables x and y, additional equations
are obtained involving the derivatives
∂z ∂z
p= , q= .
∂x ∂y
These equations can then be combined to eliminate the arbitrary constants, resulting in a PDE
involving x, y, z, and its partial derivatives. In general, the number of differentiations required
corresponds to the number of constants that must be removed.
A similar idea applies when the relation contains an arbitrary function, such as z = f (x + y)
or z = f (xy). In this case, the equation is differentiated with respect to the independent variables to
produce expressions involving the derivative of the arbitrary function. By eliminating the function
and its derivatives from these relations, one obtains a PDE satisfied by z(x, y). This procedure
is important because it shows how families of functions depending on parameters or arbitrary
functions can lead naturally to differential equations describing the entire family.
A useful rule of thumb for the formation of partial differential equations (PDEs) is that
a PDE can be obtained by eliminating arbitrary constants or arbitrary functions from a given
relation involving a dependent variable z = z(x, y). If the relation contains n arbitrary constants,
one generally differentiates the relation with respect to the independent variables x and y sufficiently
many times to obtain enough equations to eliminate those constants. The resulting equation, which
involves x, y, z, and its partial derivatives is a PDE whose order typically corresponds to the number
of independent differentiations required.
When the relation contains an arbitrary function, the idea is similar: the relation is partially
differentiated with respect to the independent variables until expressions involving the derivative of
the arbitrary function appear. These expressions are then combined to eliminate the function and
its derivatives. The resulting equation is a PDE satisfied by all members of the original family of
functions. In practice, the general strategy is therefore: differentiate the given relation with respect
to the independent variables, introduce the standard notation for partial derivatives, and eliminate
the unwanted constants or functions to obtain a differential equation involving only x, y, z, and its
derivatives. We shall describe this in the following series of examples.
Example 5.1.1. Suppose
z = ax + by.
Then we easily calculate p = a and q = b. Hence we see that z = px + qy is the required PDE.
Example 5.1.2. If
z = ax2 + by 2 ,
then we calculate p = 2ax and q = 2by. Hence ax2 = px/2 and by 2 = qy/2. Substituting this back
in the original equation gives us 2z = px + qy as the required PDE.
Example 5.1.3. If
z = f (x + y)
for some arbitrary function f , then we calculate p = f → (x+y) and q = f → (x+y). Hence the required
PDE is p = q.
Example 5.1.4. If
z = f (x2 + y 2 ),
then we calculate p = 2xf → (x2 + y 2 ) and q = 2yf → (x2 + y 2 ). Hence yp = xq and that is our required
PDE.

76
Complex analysis and diff. eq. CHAPTER 5. INTRODUCTION TO PDES

Example 5.1.5. If
z = f (xy),
then p = yf → (xy) and q = xf → (xy). Hence xp = yq is the required PDE.

5.2 Solution of first order PDEs

A linear partial differential equation is one in which the dependent variable z = z(x, y) and
its partial derivatives appear only to the first power and are not multiplied together. In general, a
linear PDE can be written in the form
A(x, y) zxx + B(x, y) zxy + C(x, y) zyy + D(x, y) zx + E(x, y) zy + F (x, y) z = G(x, y),
where the coefficients depend only on the independent variables. If the right–hand side satisfies
G(x, y) = 0, the equation is called homogeneous; otherwise it is nonhomogeneous. For example,
zxx + zyy = 0
is a linear homogeneous PDE (Laplace’s equation), while
zxx + zyy = x + y
is linear but nonhomogeneous.
The general solution of a PDE contains arbitrary constants or arbitrary functions and rep-
resents a family of solutions. In the case of a nonhomogeneous linear PDE, it is often useful to
express the solution as the sum of two parts: a particular integral, which is any single solution of
the nonhomogeneous equation, and the general solution of the associated homogeneous equation.
An important property of linear homogeneous PDEs is the superposition principle. This
principle states that if z1 (x, y) and z2 (x, y) are solutions of a linear homogeneous PDE, then any
linear combination
z = c 1 z1 + c 2 z2
is also a solution, where c1 and c2 are constants. This property allows more complicated solutions
to be constructed by combining simpler ones and plays a fundamental role in many analytical
methods for solving PDEs.

5.3 Standard types of first order PDEs

5.3.1 The equation F (p, q) = 0
Some simple first order PDEs can be solved easily. Suppose that the PDE does not involve the
dependent variable z or the independent variables x, y. That is suppose that the PDE is given in
the form F (p, q) = 0. To solve such a PDE, we guess the solution of the form z = ax + by + c.
It follows that p = a and q = b, so that F (a, b) = 0. Rewriting b in terms of a, say b = g(a), we
obtain the general solution as
z = ax + g(a)y + c.
Example 5.3.1. Suppose we want to solve
p + q = 1.
Assume z = ax + by + c. Then p = a, q = b. Substitute into the PDE a + b = 1 and hence b = 1 − a,
giving us the solution
z = ax + (1 − a)y + c.

77
Complex analysis and diff. eq. CHAPTER 5. INTRODUCTION TO PDES

Example 5.3.2. Suppose we want to solve p2 + q 2 = 1. Assume z = ax + by + c. Then p = a, q = b.

Substituting gives a2 + b2 = 1. Let a = cos θ, b = sin θ. Hence the solution is
z = x cos θ + y sin θ + c.

5.3.2 Clairut’s equation

A partial differential equation of the form z = px + qy + f (p, q) is called the Clairut’s equation.
A solution is of the form
ax + by + f (a, b)
where a, b are arbitrary constants.

5.3.3 The variable separable equation

A PDE that can be written in the form f (x, p) = g(y, q) is called a variable separable equation.
It can be solved as follows. Suppose that f (x, p) = a and g(y, q) = a. Rearranging the equation we
can write p = φ(x, a) and q = ψ(y, a). Integrating with respect x and y respectively, we see that
. .
z = φ(x, a)dx = ψ(y, a)dy.

The first integral involves an arbitrary constant that depends only on y. Similarly the second
integral involves an arbitrary constant that depends only on x. Equating them together and
eliminating arbitrary constants gives them the complete solution.
Example 5.3.3. Suppose we want to solve the differential equation p − x2 = q + y 2 . This is a variable
separable equation. We get that p = x2 + a and q = a − y 2 . Hence,
x3 y3
z= + ax + b1 (y) = ay − + b2 (x).
3 3
Rearranging, we see that the general solution is
1
z = a(x + y) + (x3 − y 3 ) + b.
3

5.3.4 The equation F (z, p, q) = 0

Suppose the original PDE involves z, p, q but not the independent variables, that is of the form
F (z, p, q) = 0.
Such a PDE can be solved by converting it into an ODE using the substitution z = z(u) where
u = x + ay. Thus p = z → (u) and q = az → (u). Thus the original PDE becomes an ODE of the form
( )
dz dz
F z, , a = 0.
du du
This ODE can be solved using techniques from the previous chapters and the solution in terms of
x, y can be obtained by substituting back for u.
Example 5.3.4. Suppose z = p2 + q 2 . Converting this back into an ODE in terms of u = x + ay, we
/ dz 02
get z = (1 + a2 ) du . Rewriting the equation, we get
2
dz z
= .
du 1 + a2
9
Solving this ODE gives us 2 (1 + a2 )z = u + b or 4(1 + a2 )z = (x + ay + b)2 .

78
Complex analysis and diff. eq. CHAPTER 5. INTRODUCTION TO PDES

5.4 Lagrange equation and Lagrange’s method

A first-order linear PDE of Lagrange type is given by

P (x, y, z) p + Q(x, y, z) q = R(x, y, z),

where p = ∂z
∂x and q = ∂y .
∂z
To solve, we form the auxiliary equations

dx dy dz
= = .
P Q R
These arise by introducing an auxiliary parameter s, so that

dx dy dz
= P, = Q, = R,
ds ds ds
which represent characteristic curves along which the PDE reduces to ODEs. Solving the auxiliary
equations yields two independent integrals

u(x, y, z) = c1 , v(x, y, z) = c2 .

Hence, the general solution is

F (u, v) = 0,
where F is an arbitrary function.
Remark 5.4.1. The parameter s is not part of the final solution; it is only used to parametrize
characteristic curves and derive the auxiliary system.
Remark 5.4.2. It is not always possible to obtain an explicit solution expressing z as a function of
x and y. Instead, solutions of a first-order PDE typically represent a family of surfaces in three-
dimensional space. From the auxiliary equations, we obtain two independent integrals u(x, y, z) =
c1 and v(x, y, z) = c2 . If we allow the constants c1 and c2 to vary independently, the resulting
equations describe a two-parameter family of curves, which is not the desired outcome. To obtain
a family of surfaces (i.e., a one-parameter family), we combine these integrals into a single relation
F (u, v) = 0, where F is an arbitrary function. This functional dependence restricts c1 and c2 ,
ensuring the solution represents a family of surfaces rather than curves.

5.5 Additional Reading : Geometric idea behind Lagrange’s method

and the Charpit’s method
A useful way to interpret solutions of partial differential equations is through geometry. Consider
a function u = u(x, y) of two independent variables. The graph of u is a surface in R3 , given by

{(x, y, u(x, y)) | (x, y) ∈ Ω ⊂ R2 }.

Thus, solving a PDE for u amounts to finding a surface in R3 whose geometry satisfies certain
prescribed conditions.
More generally, if we consider two dependent variables u = u(x, y) and v = v(x, y), then the
solution defines a parametrized surface in R4 :

(x, y) +→ (x, y, u(x, y), v(x, y)).

79
Complex analysis and diff. eq. CHAPTER 5. INTRODUCTION TO PDES

In this setting, a system of PDEs imposes constraints on the tangent planes of this surface, i.e., on
the partial derivatives ux , uy , vx , vy . Geometrically, the PDE specifies how the surface is allowed to
tilt or bend in the ambient space.
For first-order PDEs, one can interpret the equation as prescribing a field of admissible tangent
planes. A solution corresponds to a surface whose tangent plane at each point lies in this prescribed
family. For higher-order PDEs, the constraints involve curvature and higher-order geometric prop-
erties of the surface. In summary, solutions of PDEs can be viewed as geometric objects (surfaces or
submanifolds), and the differential equation encodes conditions on their local geometric structure.
From the geometric perspective above, a first-order PDE can be viewed as prescribing a family
of admissible tangent planes to a surface. A particularly important class is given by Lagrange’s
linear equation, which has the form

a(x, y, u) ux + b(x, y, u) uy = c(x, y, u).

Geometrically, this equation specifies a direction field in the extended space (x, y, u): at each point,
the vector
(a(x, y, u), b(x, y, u), c(x, y, u))
determines a direction along which the solution surface must be tangent.
A solution surface can therefore be constructed as a union of curves whose tangent vectors agree
with this field. These curves are called characteristic curves, and they are obtained by solving the
system of ordinary differential equations
dx dy du
= a(x, y, u), = b(x, y, u), = c(x, y, u).
ds ds ds
Equivalently, this system is often written in symmetric form as
dx dy du
= = .
a(x, y, u) b(x, y, u) c(x, y, u)
The key idea is that along each characteristic curve, the PDE reduces to an ordinary differential
equation. By integrating this system, one obtains two independent first integrals

φ(x, y, u) = C1 , ψ(x, y, u) = C2 ,

which represent invariant quantities along characteristics.

The general solution of Lagrange’s equation is then given implicitly by an arbitrary relation
between these integrals: / 0
F φ(x, y, u), ψ(x, y, u) = 0,
where F is an arbitrary function. Geometrically, this describes a family of surfaces formed by
combining characteristic curves in all possible ways consistent with the given direction field.
Lagrange’s linear equation represents a special class of first-order PDEs in which the derivatives
ux and uy appear linearly. In geometric terms, this corresponds to a situation where the admissible
tangent directions depend only on the point (x, y, u) and not on the orientation of the tangent
plane itself.
However, many important first-order PDEs arise in a more general nonlinear form,

F (x, y, u, p, q) = 0, where p = ux , q = uy .

In this case, the equation no longer prescribes a direction field in (x, y, u)-space alone. Instead,
it defines a constraint in the jet space (x, y, u, p, q), relating both the position and the first-order

80
Complex analysis and diff. eq. CHAPTER 5. INTRODUCTION TO PDES

derivatives. Geometrically, this means that the PDE specifies a family of admissible tangent planes
rather than just tangent directions.
To construct a solution surface, one must therefore determine not only how (x, y, u) evolve, but
also how the slopes (p, q) vary consistently along the surface. This leads to an extended system of
differential equations governing curves in the higher-dimensional space (x, y, u, p, q).
The method that accomplishes this is known as Charpit’s method. It generalizes the idea of
characteristics by introducing a system of ordinary differential equations that simultaneously tracks
the evolution of the variables (x, y, u, p, q). By integrating this system, one obtains the so-called
complete integral of the PDE, from which the general solution can be derived.
Thus, Charpit’s method can be viewed as a natural extension of the geometric approach used
for Lagrange’s equation, adapted to handle fully nonlinear first-order PDEs by lifting the problem
to a higher-dimensional geometric setting.
Consider a general first-order nonlinear PDE
F (x, y, u, p, q) = 0, where p = ux , q = uy .
Unlike the linear case, the equation now constrains both the position (x, y, u) and the slopes (p, q),
so we work in the extended space (x, y, u, p, q).
Charpit’s method constructs solution surfaces by integrating a system of ordinary differential
equations (the Charpit equations) given by
dx dy du
= Fp , = Fq , = pFp + qFq ,
ds ds ds
dp dq
= −Fx − pFu , = −Fy − qFu ,
ds ds
where subscripts denote partial derivatives of F .
Equivalently, this system is often written in symmetric form:
dx dy du dp dq
= = = = .
Fp Fq pFp + qFq −Fx − pFu −Fy − qFu
To solve the PDE using Charpit’s method, one proceeds as follows:
1. Form the function F (x, y, u, p, q) = 0 and compute the partial derivatives Fx , Fy , Fu , Fp , Fq .
2. Write down the Charpit auxiliary equations:
dx dy du dp dq
= = = = .
Fp Fq pFp + qFq −Fx − pFu −Fy − qFu
3. Solve the system to obtain relations (first integrals) between x, y, u, p, q. Typically, one first
finds relations involving p and q.
4. Use these relations together with the original PDE to eliminate p and q, expressing them in
terms of x, y, u and arbitrary constants.
5. Substitute p = ux and q = uy to obtain a reduced system of first-order equations.
6. Integrate this system to obtain a family of solutions depending on two arbitrary constants.
This yields the complete integral
u = u(x, y; a, b).
7. If required, eliminate the constants a, b from
u = u(x, y; a, b), φ(a, b) = 0,
to obtain the general solution.

81