7More partial differential equations
IB Methods
7.4 Green’s functions for PDEs on R
n
Green’s functions for the heat equation
Suppose φ : R
n
× [0, ∞) → R solves the heat equation
∂
t
φ = D∇
2
φ,
where D is the diffusion constant, subject to
φ|
R
n
×{0}
= f.
To solve this, we take the Fourier transform of the heat equation in the spatial
variables. For simplicity of notation, we take
n
= 1. Writing
F
[
φ
(
x, t
)] =
˜
φ
(
k, t
),
we get
∂
t
˜
φ(k, t) = −Dk
2
˜
φ(k, t)
with the boundary conditions
˜
φ(k, 0) =
˜
f(k).
Note that this is just an ordinary differential equation in
t
, and
k
is a fixed
parameter for each k. So we have
˜
φ(k, t) =
˜
f(k)e
−Dk
2
t
.
We now take the inverse Fourier transform and get
φ(x, t) =
1
2π
Z
R
e
ikx
h
˜
f(k)e
−Dk
2
t
i
dk
This is the inverse Fourier transform of a product. This thus gives the convolution
φ(x, t) = f ∗ F
−1
[e
−Dk
2
t
].
So we are done if we can find the inverse Fourier transform of the Gaussian.
This is an easy exercise on example sheet 3, where we find
F[e
−a
2
x
2
] =
√
π
a
e
−k
2
/4a
2
.
So setting
a
2
=
1
4Dt
,
So we get
F
−1
[e
−Dk
2
t
] =
1
√
4πDt
exp
−
x
2
4Dt
We shall call this
S
1
(
x, t
), where the subscript 1 tells us we are in 1+1 dimensions.
This is known as the fundamental solution of the heat equation. We then get
φ(x, t) =
Z
∞
−∞
f(y)S
1
(x − y, t) dy
Example. Suppose our initial data is
f(x) = φ
0
δ(x).
So we start with a really cold room, with a huge spike in temperature at the
middle of the room. Then we get
φ(x, t) =
φ
0
√
4πDt
exp
−
x
2
4Dt
.
What this shows, as we’ve seen before, is that if we start with a delta function,
then as time evolves, we get a Gaussian that gets shorter and broader.
Now note that if we start with a delta function, at
t
= 0, everywhere outside
the origin is non-zero. However, after any infinitely small time
t
,
φ
becomes
non-zero everywhere, instantly. Unlike the wave equation, information travels
instantly to all of space, i.e. heat propagates arbitrarily fast according to this
equation (of course in reality, it doesn’t). This fundamentally, is because the
heat equation is parabolic, and only has one characteristic.
Now suppose instead φ satisfies the inhomogeneous, forced, heat equation
∂
t
φ − D∇
2
φ = F (x, t),
but with homogeneous initial conditions
φ|
t=0
= 0. Physically, this represents
having an external heat source somewhere, starting at zero.
Note that if we can solve this, then we have completely solved the heat
equation. If we have an inhomogeneous equation and an inhomogeneous initial
condition, then we can solve this forced problem with homogeneous boundary
conditions to get
φ
F
; and solve the unforced equation with homogeneous equation
to get φ
H
. Then the sum φ = φ
F
+ φ
H
solves the full equation.
As before, we take the Fourier transform of the forced equation with respect
to the spacial variables. As before, we will just do it in the cases with one spacial
dimension. We find
∂
t
˜
φ(k, t) + Dk
2
˜
φ(k, t) =
˜
F (k, t),
with the initial condition
˜
φ(k, 0) = 0.
As before, we have reduced this to a first-order ordinary differential equation in
t. Using an integrating factor, we can rewrite this as
∂
∂t
[e
Dk
2
t
˜
φ(k, t)] = e
Dk
2
t
˜
F (k, t).
The solution is them
˜
φ(k, t) = e
−Dk
2
t
Z
t
0
e
Dk
2
u
˜
F (k, u) du.
We define the Green’s function G(x, t; y, τ)) to be the solution to
[∂
t
− D∇
2
x
]G(x, t; y, τ) = δ(x − y)δ(t − τ ).
So the Fourier transform with respect to x gives
˜
G(k, t, y, τ) = e
−Dk
2
t
Z
t
0
e
Dk
2
u
e
iky
δ(t − τ) du,
where e
iky
is just the Fourier transform of δ(x − y). This is equal to
˜
G(k, t; y, τ) =
(
0 t < τ
e
−iky
e
−Dk
2
(t−τ)
t > τ
= Θ(t − τ)e
−iky
e
−Dk
2
(t−τ)
.
Reverting the Fourier transform, we get
G(x, t; y, τ) =
Θ(t − τ)
2π
Z
R
e
ik(x−y)
e
−Dk
2
(t−τ)
dx
This integral is just the inverse Fourier transform of the Gaussian with a phase
shift. So we end up with
G(x, t; y, τ) =
Θ(t − τ)
p
4πD(t − τ)
exp
−
(x − y)
2
4D(t − τ)
= Θ(t − τ)S
1
(x − y; t − τ).
The solution we seek is then
φ(x, t) =
Z
t
0
Z
R
F (y, τ)G(x, t; y, τ ) dy dτ.
It is interesting that the solution to the forced equation involves the same function
S
1
(
x, t
) as the homogeneous equation with inhomogeneous boundary conditions.
In general, S
n
(x, t) solves
∂S
n
∂t
− D∇
2
S
n
= 0
with boundary conditions S
n
(x, 0) = δ
(n)
(x − y), and we can find
S
n
(x, t) =
1
(4πDt)
n/2
exp
−
|x − y|
2
4Dt
.
Then in general, given an initial condition φ|
t=0
= f(x), the solution is
φ(x, t) =
Z
f(y)S(x − y, t) d
n
y.
Similarly, G
n
(x, t; y, t) solves
∂G
n
∂t
− D∇
2
G
n
= δ(t − τ)δ
(n)
(x − y),
with the boundary condition G
n
(x, 0; y, τ ) = 0. The solution is
G(x, t; y, τ ) = Θ(t − τ )S
n
(x − y, t − τ ).
Given our Green’s function, the general solution to the forced heat equation
∂φ
∂t
− D∇
2
φ = F (x, t), φ(x, 0) = 0
is just
φ(x, t) =
Z
∞
0
Z
R
n
F (y, τ )G(x, t; y, τ) d
n
y dτ
=
Z
t
0
Z
R
n
F (y, τ )G(x, t; y, τ) d
n
y dτ.
Duhamel noticed that we can write this as
φ(x, t) =
Z
t
0
φ
F
(x, t; τ) dτ,
where
φ
F
=
Z
R
F (y, t)S
n
(x − y, t − τ ) d
n
y
solves the homogeneous heat equation with φ
F
|
t=τ
= F (x, τ).
Hence in general, we can think of the forcing term as providing a whole
sequence of “initial conditions” for all
t >
0. We then integrate over the times
at which these conditions were imposed to find the full solution to the forced
problem. This interpretation is called Duhamel’s principle.
Green’s functions for the wave equation
Suppose φ : R
n
× [0, ∞) → C solves the inhomogeneous wave equation
∂
2
φ
∂t
2
− c
2
∇
2
φ = F (x, t)
with
φ(x, 0) =
∂
∂t
φ(x, 0) = 0.
We look for a Green’s function G
n
(x, t; y, τ ) that solves
∂G
n
∂t
− c
2
∇
2
G
n
= δ(t − τ)δ
(n)
(x − y) (∗)
with the same initial conditions
G
n
(x, 0, y, τ ) =
∂
∂t
G
n
(x, 0, y, τ ) = 0.
Just as before, we take the Fourier transform of this equation with respect the
spacial variables x. We get
∂
∂t
˜
G
n
+ c
2
|k|
2
˜
G
n
= δ(t − τ)e
−ik·y
.
where
˜
G
n
=
˜
G
n
(k, t, y, τ ).
This is just an ordinary differential equation from the point of view of
t
, and
is of the same type of initial value problem that we studied earlier, and the
solution is
˜
G
n
(k, t, y, τ ) = Θ(t − τ )e
−ik·y
sin |k|c(t − τ)
|k|c
.
To recover the Green’s function itself, we have to compute the inverse Fourier
transform, and find
G
n
(x, t; y, τ ) =
1
(2π)
n
Z
R
n
e
ik·x
Θ(t − τ)e
−ik·y
sin |k|c(t − τ)
|k|c
d
n
k.
Unlike the case of the heat equation, the form of the answer we get here does
depend on the number of spatial dimensions
n
. For definiteness, we look at the
case where
n
= 3, since our world (probably) has three dimensions. Then our
Green’s function is
G(x, t; y, τ ) =
Θ(t − τ)
(2π)
3
c
Z
R
3
e
ik·(x−y)
sin |k|c(t − τ)
|k|
d
3
k.
We use spherical polar coordinates with the
z
-axis in
k
-space aligned along the
direction of x − y. Hence k · (x − y) = kr cos θ, where r = |x − y| and k = |k|.
Note that nothing in our integral depends on
ϕ
, so we can pull out a factor
of 2π, and get
G(x, t; y, τ ) =
Θ(t − τ)
(2π)
2
c
Z
∞
0
Z
π
0
e
ikr cos θ
sin kc(t − τ)
k
k
2
sin θ dθ dk
The next integral to do is the
θ
integral, which is straightforward since it is an
exact differential. Setting α = c(t − τ), we get
=
Θ(t − τ)
(2π)
2
c
Z
∞
0
e
ikr
− e
−ikr
ikr
sin kc(t − τ)
k
k
2
dk
=
Θ(t − τ)
(2π)
2
icr
Z
∞
0
e
ikr
sin kα dk −
Z
∞
0
e
−ikr
sin kα dk
=
Θ(t − τ)
2πicr
1
2π
Z
∞
−∞
e
ikr
sin kα dk
=
Θ(t − τ)
2πicr
F
−1
[sin kα],
Now recall F[δ(x − α)] = e
−ikα
. So
F
−1
[sin kα] = F
−1
e
ikα
− e
−ikα
2i
=
1
2i
[δ(x + α) − δ(x − α)].
Hence our Green’s function is
G(x, t; y, τ ) = −
Θ(t − τ)
4πc|x − y|
h
δ
|x − y| + c(t − τ )
− δ
|x − y| − c(t − τ )
i
.
Now we look at our delta functions. The step function is non-zero only if
t > τ
.
Hence
|x − y|
+
c
(
t − τ
) is always positive. So
δ
(
|x − y|
+
c
(
t − τ
)) does not
contribute. On the other hand,
δ
(
|x −y|−c
(
t −τ
)) is non-zero only if
t > τ
. So
Θ(
t − τ
) is always positive in this region. So we can write our Green’s function
as
G(x, t; y, τ ) =
1
4πc
1
|x − y|
δ(|x − y | − c(t − τ)).
As always, given our Green’s function, the general solution to the forced equation
∂
2
φ
∂t
2
− c
2
∇
2
φ = F (x, t)
is
φ(x, t) =
Z
∞
0
Z
R
3
F (y, τ )
4πc|x − y|
δ(|x − y | − c(t − τ)) d
3
y dτ.
We can use the delta function to do one of the integrals. It is up to us which
integral we do, but we pick the time integral to do. Then we get
φ(x, t) =
1
4πc
2
Z
R
3
F (y, t
ret
)
|x − y|
d
3
y,
where
t
ret
= t −
|x − y|
c
.
This shows that the effect of the forcing term at some point
y ∈ R
3
affects
the solution
φ
at some other point
x
not instantaneously, but only after time
|x − y|/c
has elapsed. This is just as we saw for characteristics. This, again,
tells us that information travels at speed c.
Also, we see the effect of the forcing term gets weaker and weaker as we
move further away from the source. This dispersion depends on the number of
dimensions of the space. As we spread out in a three-dimensional space, the
“energy” from the forcing term has to be distributed over a larger sphere, and
the effect diminishes. On the contrary, in one-dimensional space, there is no
spreading out to do, and we don’t have this reduction. In fact, in one dimensions,
we get
φ(x, t) =
Z
t
0
Z
R
F (y, τ)
Θ(c(t − τ) − |x − y|)
2c
dy dτ.
We see there is now no suppression factor at the bottom, as expected.