1. business and industrial

Green’s functions, formulas and representations
Suppose that we want to solve a linear, inhomogeneous equation of the
form
Lu(x) = f (x) + homogeneous boundary conditions.
(1)
Here u, f are functions whose domain is Ω, which could be either a bounded
or unbounded domain in any dimension. The boundary conditions we will
consider will be Dirichlet, Neumann, or some mixture.
In this discussion, L will be self-adjoint with respect to the usual L2
inner product using the imposed homogeneous boundary conditions. In
addition, we will assume that equation (1) admits a unique solution for any
(sufficiently smooth) right hand side f . Both assumptions can be relaxed,
by the way, using adjoint- or modified- Green’s functions, although these are
not pursued here.
1
The delta function and distributions
There is a great need in differential equations to define objects that arise
as limits of functions and behave like functions under integration but are
not, properly speaking, functions themselves. These objects are sometimes
called generalized functions or distributions. The most basic one of these is
the so-called δ-function.
For each > 0, define the family of ordinary functions
1
2 2
δ (x) = √ e−x / .
π
(2)
When is small, the graph of δ is essentially just a spike at x = 0, but the
area under the curve is exactly one. Thus for any continuous and bounded
f (x),
Z ∞
Z ∞
lim
f (x)δ (x − x0 )dx = lim f (x0 )
δ (x − x0 )dx = f (x0 ).
→0 −∞
→0
−∞
On the other hand, taking the limit → 0 inside the integral makes no sense:
the limit of δ is not a function at all! To get around this, we define the object
δ to act as follows:
Z ∞
f (x)δ(x − x0 )dx = f (x0 ).
(3)
−∞
1
Informally speaking, the δ-function “picks out” the value of the function
f (x) at one point.
There are δ functions for higher dimensions also. We define the ndimensional δ function as
Z
f (x)δ(x − x0 )dx = f (x0 ),
Rn
where x, x0 ∈ Rn . Sometimes we write this higher dimensional delta function as a product of one dimensional ones δ(x) = δ(x) · δ(y) · . . ..
1.1
A more precise definition
To be concrete about distributions, one needs to talk about how they “act”
on smooth functions. Note that (3) can be thought of as a linear mapping
from smooth functions f (x) to real numbers by the process of evaluating
them at x0 . Linear mappings from a vector space to the real numbers are
often called linear functionals.
Now we come to the precise definition. A distribution is a continuous linear functional on the set of infinitely differentiable functions with
bounded support (this vector space of functions is denoted C0∞ or simply
D). We can write d[φ] : D → R to represent such a map: for any input
function φ, d[φ] gives us a number.
The class of distributions includes just plain old functions in the following sense. For some integrable function g(x), we define the corresponding
linear functional to be
Z ∞
g[φ] =
g(x)φ(x)dx.
−∞
(notice the subtlety of notation: g(x) means the function evaluated at x,
whereas the bracket g[φ] means multiply by function φ and integrate). In
fact, all distributions can be approximated by smooth, ordinary functions.
This means that, for any distribution d, there exists a sequence of smooth
functions dn (x) ∈ D so that
Z
dn (x)φ(x)dx, for all φ ∈ D.
d[φ] = lim
n→∞
For example, the sequence δ that we first investigated are smooth approximations to the delta function distribution.
2
We can now define what it means to integrate a distribution (notated as
if it were a regular function d(x)) simply by
Z ∞
d(x)φ(x)dx = d[φ], for any φ ∈ D.
−∞
This is consistent with the formula (3) since δ(x) maps a function φ onto its
value at zero.
Here are a couple examples. A linear combination of two delta functions such as d = 3δ(x−1)+2δ(x) defines a distribution. The corresponding
linear functional is
Z ∞
d(x)φ(x)dx.
d[φ] = 3φ(1) + 2φ(0) =
−∞
The linear operation
d[φ] = φ0 (0)
is a linear map (check it!), taking a smooth function and returning the value
of its derivative at zero. A sequence of smooth functions that approximates
this is −δ0 (x), since integration by parts gives
Z
Z
0
lim −δ (x)φ(x)dx = lim δ (x)φ0 (x)dx = φ0 (0),
→0
→0
1.2
Distributions as derivatives
One useful aspect of distributions is that they make sense of derivatives of
functions which are non-smooth or even unbounded. Suppose that g(x)
is an integrable function, but cannot be differentiated everywhere in its
domain. It can make sense, however, to talk about integrals involving g 0 .
Though integration by parts doesn’t technically hold in the usual sense, for
φ ∈ D we can define
Z ∞
Z ∞
g 0 (x)φ(x)dx ≡ −
g(x)φ0 (x)dx.
−∞
−∞
Notice that the expression on the right makes perfect sense as a usual integral. We define the distributional derivative of g(x) to be the distribution
g 0 [φ] so that
Z ∞
g 0 [φ] ≡ −
g(x)φ0 (x)dx, for all φ ∈ D.
−∞
3
Notice that g does not even need to be a function for this definition to make
sense. If g is merely a distribution, then its derivative is the rule given by
g 0 [φ] ≡ −g[φ0 ].
For example, if H(x) is the step function which is zero when x < 0 and
one when x > 0, the distributional derivative H 0 [φ] of H is the rule
Z ∞
Z ∞
Z ∞
δ(x)φ(x)dx.
φ0 (x)dx = φ(0) =
H(x)φ0 (x)dx = −
H 0 [φ] = −
−∞
−∞
0
Therefore the delta function is the distributional derivative of the unit step
function.
1.3
The distributional Laplacian
In higher dimensions, one can make similar definitions of distributional
derivatives by using Green’s identities. For smooth functions u(x) and φ ∈
D, one has
Z
Z
(∆u)φdx =
u ∆φdx,
Rn
Rn
since the boundary terms in the Green identity vanish. This motivates a
definition of the distributional Laplacian for non-smooth functions, which is
a distribution with linear functional
Z
(∆u)[φ] =
u ∆φdx.
Rn
As an example,
p let’s compute the distributional Laplacian of f (x, y) =
ln(sinh r), r = x2 + y 2 . The definition tells us we want to know what rule
is implied by
Z 2π Z ∞
h1
i
1
∆f [φ] = f [∆φ] =
ln(sinh r) (rφr )r + 2 φθθ rdrdθ.
r
r
0
0
The integral of the second term is
Z ∞ Z 2π
1
φθθ dθ dr = 0,
r
0
0
R 2π
since the θ integral is zero. Now we can integrate the term inside 0 dθ by
parts twice:
Z ∞
∞ Z ∞
ln(sinh r)(rφr )r dr = r ln(sinh r)φr −
r coth rφr dr =
0
0
0
∞ Z ∞
− φ(r)r coth r +
(r coth r)r φdr.
0
0
4
By using L’Hospital’s rule, the limit φ(r)r coth r → φ(0) as r → 0. Then it
follows
Z ∞ Z 2π
Z
1
1
f [∆φ] = 2πφ(0)+
(r coth r)r φ rdrdθ = 2πφ(0)+
(r coth r)r φdr.
r
2
0
0
R r
Notice that to turn the iterated integral back into an area integral, a factor of
r was needed. Therefore the distribution Laplacian is 2πδ(x) + 1r (r coth r)r ;
it is the sum of a regular function and one which can only be regarded as a
distribution.
2
Green’s functions
Recall that differential operators often have inverses that are integral operators. So for equation (1), we might expect a solution of the form
Z
u(x) =
G(x; x0 )f (x0 )dx0
(4)
Ω
If such a representation exists, the kernel of this integral operator G(x; x0 )
is called the Green’s function.
It is useful to give a physical interpretation of (4). We think of u(x) as
the response at x to the influence given by f (x). For example, if the problem involved elasticity, u might be the displacement subject to an external
force f . If this were a heat equation, u might be the temperature subject
to a heat source described by f . The integral can be though of as the sum
(or superposition) over influences created by sources at each x0 . For this
reason, G is sometimes called the influence function.
2.1
Relationship to the delta function
Part of the problem with the definition (4) is that it doesn’t tell us how to
construct G. It is useful to imagine what happens when f is a point source
at xi ∈ Ω, i.e. f (x) = δ(x − xi ). Plugging into (4) we learn that the solution
to
Lu(x) = δ(x − xi ) + homogeneous boundary conditions
(5)
is just
Z
G(x; x0 )δ(x0 − xi )dx0 = G(x; xi ).
u(x) =
(6)
Ω
In other words, we find that the Green’s function G(x; x0 ) formally satisfies
Lx G(x; x0 ) = δ(x − x0 )
5
(7)
(the subscript means that the linear operator acts on x, not x0 ). This equation says that G(x; x0 ) is the influence felt at point x due to a source at point
x0 .
Equation (7) is a more useful way of defining G since we can in many
cases solve this “almost” homogeneous equation, either by direct integration or using Fourier techniques. In particular, (7) can be rewritten as
Lx G(x, x0 ) = 0,
when x 6= x0 .
(8)
To account for the δ-function, we can formally integrate both sizes of Lx G(x; x0 ) =
δ(x − x0 ) on any region containing x0 . It is usually sufficient to allow these
regions to be balls Br (x0 ) = {x||x − x0 | < r}, so that
Z
Lx G(x, x0 )dx = 1.
(9)
Br (x0 )
Equation (8) is a homogeneous equation with a “hole” in the domain at
x0 . Equation (9) is called the normalization condition, and it is used to get
the “size” of the singularity at x0 correct. In one dimension, this condition
takes on a slightly different form (see below).
In addition to (8-9), G must also satisfy the same type of homogeneous
boundary conditions that the solution u does in the original problem. The
reason for this is straightforward. Take, for example, the case of a homogeneous Dirichlet boundary condition u = 0 for x ∈ ∂Ω. For any point x on
the boundary, it must be the case that
Z
G(x; x0 )f (x0 )dx0 = 0.
(10)
Ω
Since this must be true for any choice of f , it follows that G(x; x0 ) = 0 for
boundary points x (note that x0 is treated as a constant in this respect, and
can be any point in the domain).
2.2
Jump conditions for ODEs
Instead of the integral condition (9), Green’s functions in one dimension
(for ordinary differential equations) will satisfy pointwise conditions for
their behavior when x → x0 . Suppose one has a n-th order linear equation
of the form
u(n) (x) + F (u(n−1) (x), u(n−2) (x), . . .) = f (x),
6
where F is some expression involving lower order derivatives. The Green’s
function G(x; x0 ) formally satisfies
G(n) + F (G(n−1) , G(n−2) , . . .) = δ(x − x0 ),
where G(n) = ∂ n /∂xn . Formal integration of both sides gives
Gn−1 = H(x − x0 ) + some continuous function
This means that something special happens at x0 : the n − 1-th derivative is
not continuous, but suffers a discontinuous jump there. Integrating again
shows that Gn−2 is continuous.
These facts can be summarized as “jump” and continuity conditions at
x0 :
lim
x→x+
0
∂ n−1 G
∂ n−1 G
−
lim
= 1,
∂xn−1 x→x−
∂xn−1
0
lim
x→x+
0
∂mG
∂mG
−
lim
= 0,
∂xm x→x−
∂xm
0
m ≤ n−2.
(11)
For first order equations, (11) means that G itself must have a jump discontinuity. For second order equations, G is continuous but its derivative
has a jump discontinuity.
2.3
Distributional interpretation for the Laplacian
Let’s specialize to the most important linear operator, the Laplacian. We
would like give precise meanings to the equations
u00 (x) = δ(x − x0 ),
∆u(x) = δ(x − x0 ),
x∈R
(12)
n
x∈R , n≥2
(13)
For the one dimensional situation (12), we can interpret u00 as a distribution acting on φ ∈ D as
Z ∞
u00 [φ] =
u(x)φ00 (x)dx = φ(x0 )
−∞
since δ(x − x0 ) is just the distribution which picks out the value of φ at x0 .
If we assume u(x) is a smooth function everywhere but x0 , integration by
parts twice is justified providing the integral is split first:
Z x0
Z ∞
00
φ(x0 ) =
u(x)φ (x)dx +
u(x)φ00 (x)dx
(14)
−∞
x0
Z x0
Z ∞
=
u00 (x)φ(x)dx +
u00 (x)φ(x)dx + φ(x0 )[u0 ]x0 − φ0 (x0 )[u]x0 ,
−∞
x0
(15)
7
where [f ]z denotes the difference of the right and left limits of f at z. If we
consider any φ(x) which is nonzero for x < x0 , we find that
Z x0
u00 (x)φ(x)dx = 0,
−∞
for any such φ. The only way this can happen is if u00 (x) = 0 for x < x0 .
Similarly u00 (x) = 0 for x > x0 , which means that (15) becomes
φ(x0 ) = φ(x0 )[u0 ]x0 − φ0 (x0 )[u]x0 .
Choosing φ(x0 ) = 1 and φ0 (x0 ) = 0 means [u0 ]x0 = 1, which is the jump
condition on the derivative of u in (11). Conversely, choosing φ(x0 ) = 0
and φ0 (x0 ) = 1 means [u]x0 = 0, which is the continuity condition in (11).
Now consider higher dimensional problem (13). This means that the
distributional Laplacian of u is
Z
u∆φ dx = δ[φ] = φ(0),
(16)
Rn
for every φ ∈ D. Suppose first we choose a point x0 6= 0, and let φ be any
function whose support is in a small neighborhood B of x0 . Then supposing u is a regular smooth function near x0 ,
Z
Z
φ∆u dx =
u∆φ = φ(0) = 0,
(17)
B
B
by using (16). Since this is true for any choice of φ with φ = 0 outside B,
this means that
∆u(x0 ) = 0, x0 6= 0.
(18)
Conversely, take φ so that φ = 1 on some neighborhood B containing the
origin. Using (16),
Z
Z
1 = φ(0) =
u∆φ dx =
u∆ φdx,
(19)
Rn
R/B
since the integrand is zero on B. We can now use Green’s identity on the
remaining integral, since u is a perfectly nice function away from the origin.
By virtue of (18), and the fact that ∇φ = 0 on ∂B, we are left with
Z
∇u · n
ˆ dx = 1.
(20)
∂B
One way of thinking of this statement is that the vector field ∇u has a unit
source at the origin. The two conditions (18),(20) will be needed in the next
section to find such a function u(x).
8
2.4
Three standard examples
We now actually find Green’s functions in one, two and three dimensions.
One dimension. Suppose u : R → R solves the (ordinary differential equation)
uxx = f, u(0) = 0 = u(L).
(21)
The corresponding Green’s function will solve
Gxx (x; x0 ) = 0 for x 6= x0 ,
G(0, x0 ) = 0 = G(L, x0 ),
(22)
The jump conditions are
lim Gx (x; x0 ) − lim Gx (x; x0 ) = 1, lim G(x; x0 ) − lim G(x; x0 ) = 0.
x→x+
0
x→x−
0
x→x+
0
x→x−
0
(23)
The problem in (22) actually represents two boundary value problems,
one to the left of x0 and one to the right. Their solution is
(
c1 x
x < x0
G(x; x0 ) =
(24)
c2 (x − L) x > x0
Now imposing conditions (23) gives
c1 x0 = c2 (x0 − L),
c2 − c1 = 1,
so that c1 = (x0 − L)/L and c2 = x0 /L. It follows that the solution to (21)
can be written using G as
Z x
Z L
1
u(x) =
x0 (x − L)f (x0 )dx0 +
x(x0 − L)f (x0 )dx0 .
L
0
x
Three dimensions. Now suppose u : R3 → R solves
∆u = f,
lim u(x) = 0.
|x|→∞
(25)
In this case the homogeneous “boundary” condition is actually a far-field
condition. The corresponding Green’s function therefore must solve (8),
∆x G(x; x0 ) = 0,
x 6= x0 ,
lim G(x, x0 ) = 0.
|x|→∞
The normalization condition (9) is the same as (20),
Z
∇x G(x; x0 ) · n
ˆ dx = 1.
∂B
9
(26)
(27)
for any ball B centered on x0 . Let us observe the following: if we rotate
the Green’s function about x0 , it still will solve (26-27) since the Laplace
operator is invariant under rotation. Thus G only depends on the distance
between x and x0 . We can write G = G(r), r = |x − x0 |, where in spherical
coordinates (26) is
1 2 0
(r G (r))0 = 0 if r 6= 0,
r2
lim G(r) = 0.
r→∞
(28)
This is easily integrated twice to give the general solution
G=−
c1
+ c2 ,
r
(29)
where c2 = 0 by using the far-field condition in (28). The normalization
condition (27) determines c1 . Letting B be the unit ball centered at x0
(whose surface area is 4π),
Z
Z
c1
1=
∇x G(x; x0 ) · n
ˆ dx =
dx = 4πc1 ,
2
r
∂B
∂B
so that c1 = 1/4π. Thus the Green’s function is G(x; x0 ) = −1/(4π|x − x0 |),
and the solution to (25) is
Z
f (x0 )
u(x) = −
dx30 .
R3 4π|x − x0 |
Two dimensions. In this case we want to solve
∆u = f, lim u(r, θ) − ur (r, θ)r ln r = 0.
r→∞
(30)
The far field condition looks very strange at first glance. The reason for this
is that solutions to ∆u = f will in general NOT decay to zero at large |x|,
but will grow logarithmically. Again we look for a Green’s function of the
form G = G(|x − x0 |) = G(r), so that
1
(31)
(rG0 (r))0 = 0 if r 6= 0, lim G − Gr r ln r = 0.
r→∞
r
The general solution is
G = c1 ln r + c2 ,
(32)
where c2 = 0 by using the far-field condition in (31). The normalization
condition (27) gives (where B is the unit disk)
Z
Z
1=
∇x G(x; x0 ) · n
ˆ dx =
c1 dx = 2πc1 ,
∂B
∂B
10
so that c1 = 1/2π. Thus the Green’s function is G(x; x0 ) = ln |x − x0 |/2π,
and the solution to (30) is
Z
ln |x − x0 |f (x0 ) 2
u(x) =
dx0 .
2π
2
R
It is sometimes useful to write G in polar coordinates. Using the law of
cosines for the distance |x − x0 |, one gets
G(r, θ; ro , θ0 ) =
3
1
ln(r2 + r02 − 2rr0 cos(θ − θ0 )).
4π
(33)
Boundaries and the method of images
The examples in the previous section are free space Green’s functions, since
there are no domain boundaries. Recall that the Green’s function must
satisfy all the same homogeneous boundary conditions as the PDE to be
solved for u. Unfortunately, the free space Green’s functions only satisfy a
condition at infinity. But with a little cleverness, we can still employ them
as a starting point even when the domain is bounded or semi-infinite.
3.1
Arbitrary bounded domains
Suppose that we wish to solve the Poisson equation for u : Ω → R
∆u = f (x, y),
u = 0 on ∂Ω
(34)
where Ω is a bounded, open set in R2 . We need a Green’s function G(x, x0 )
which satisfies
∆x G = δ(x − x0 ),
G(x, x0 ) = 0 when x ∈ ∂Ω.
(35)
It’s tempting to use the free space Green’s function G2 (x; x0 ) = ln |x −
x0 |/2π, which does indeed satisfy the equation in (35). Unfortunately G2 is
not zero on the boundary ∂Ω!
We should be thinking of (35) as an inhomogeneous equation, and use
the method of particular solutions. In fact, G2 is a particular solution, so if
we write G(x, x0 ) = G2 (x; x0 ) + GR (x, x0 ), then GR solves a homogeneous
equation for each x0 ∈ Ω
∆x GR = 0,
GR (x, x0 ) = −G2 (x, x0 ) when x ∈ ∂Ω.
11
(36)
The function GR is called the regular part of the Green’s function, and it
just a nice, usual solution of Laplace’s equation, without any nasty singular behavior. In other words, the desired Green’s function has exactly the
same logarithmic singularity as the free space version, but is quantitatively
different far away from x0 .
3.2
Method of images
We can often take advantage of the symmetry of domains to construct
Green’s functions using free space Green’s functions as building blocks.
This comes from a fairly simple observation about even and odd functions:
for continuously differentiable f (x) : R → R,
g(x) = f (x) − f (−x) is an odd function and g(0) = 0,
0
h(x) = f (x) + f (−x) is an even function and h (0) = 0.
(37)
(38)
Therefore subtracting mirror images of a function will satisfy a Dirichlettype boundary condition at x = 0, whereas their sum satisfies a Neumanntype boundary condition.
As an example, consider finding the Green’s function for the upper-half
space problem for u : R3 ∩ {z > 0} → R:
∆u = f,
lim u(x) = 0,
|x|→∞
(39)
u(x, y, 0) = 0.
The free space Green’s function (in Cartesian coordinates) is
G3 (x, y, z; x0 , y0 , z0 ) = −
1
p
2
4π (x − x0 ) + (y − y0 )2 + (z − z0 )2
is not zero where z = 0. On the other hand, using (37), subtracting G from
its mirror image about z = 0 should do the trick:
G(x, y, z; x0 , y0 , z0 ) = G3 (x, y, z; x0 , y0 , z0 ) − G3 (x, y, −z; x0 , y0 , z0 )
1
=
4π
1
1
p
−p
(x − x0 )2 + (y − y0 )2 + (z − z0 )2
(x − x0 )2 + (y − y0 )2 + (z + z0 )2
Let’s check that this is what we want. It is easy to see that lim|x|→∞ G = 0
and also that G(x, y, 0; x0 , y0 , z0 ) = 0. The equation formally satisfied by G
is
∆x G(x, x0 ) = δ(x − x0 ) − δ(x − x∗0 ),
(40)
12
!
.
where x∗0 = (x0 , y0 , −z0 ) is called the image source. The presence of the
extra δ-function on the right hand side in (40) is not a problem, because a
differential equation only needs to be satisfied in the domain where it is
posed. Therefore
∆x G(x, x0 ) = δ(x − x0 ),
for all x, x0 ∈ R3 ∩ {z > 0} .
A more clever use of the method of images is where Ω is a disk of radius
a, and we want (in polar coordinates) G(r, θ; r0 , θ0 ) to solve ∆G = δ(x −
x0 ) with boundary condition G(a, θ; r0 , θ0 ) = 0. The idea is to subtract
the two-dimensional free space Green’s function G2 from some carefully
chosen “reflection” across the boundary of the circle. This results (after
some guesswork!) in
2
1
r2 + r02 − 2rr0 cos(θ − θ0 )
a
G(r, θ; r0 , θ0 ) =
(41)
ln
4π
r02 r2 + a4 /r02 − 2ra2 /r0 cos(θ − θ0 )
1
= G2 (x, x0 ) − G2 (x, x∗0 ) +
ln(a/r0 ),
(42)
2π
where the image source is outside the disk at x∗0 = a2 x0 /r02 . Does this work?
First, we have ∆G = δ(x − x0 ) − δ(x − x∗0 ) as in the previous example,
which is just ∆G = δ(x − x0 ) when restricted to the disk. Evaluating G on
the boundary,
2
1
a
a2 + r02 − 2ar0 cos(θ − θ0 )
G(a, θ; ro , θ0 ) =
ln
4π
r2 a2 + a4 /r02 − 2a3 /r0 cos(θ − θ0 )
02
1
a + r02 − 2ar0 cos(θ − θ0 )
= 0.
=
ln 2
4π
r0 + a2 − 2ar0 cos(θ − θ0 )
4
Inhomogeneous boundary conditions
Remarkably, the Green’s function can be used for problems with inhomogeneous boundary conditions even though the Green’s function itself satisfies
homogeneous boundary conditions. To obtain a representation formula for the
solution, we will need a ”Green’s formula” particular to the linear operator
in question.
4.1
Green’s formulas
If a differential operator L is self-adjoint with respect to the usual L2 inner
product, then for all functions u, v satisfying the homogeneous boundary
13
conditions in problem (1),
Z
Z
(Lv)u dx − (Lu)v dx = 0.
Ω
(43)
Ω
What if u, v don’t necessarily satisfy homogeneous boundary conditions?
Then something like (43) would still be true, but terms involving boundary
values of u, v would appear:
Z
Z
(44)
(Lv)u dx − (Lu)v dx = boundary terms involving u and v.
Ω
Ω
What this formula actually looks like depends on the linear operator in
question and is known as the Green’s formula for the linear operator L.
4.2
Symmetry (reciprocity) of the Green’s function
Occasionally we need to utilize a special property of Green’s functions of
self-adjoint operators often called reciprocity. This states that G(x; x0 ) =
G(x0 ; x), which means that the associated integral operator is self-adjoint.
To show this, let v(x) = G(x; x1 ) and u(x) = G(x; x2 ). Then because of (5),
we have Lv = δ(x − x1 ) and Lu = δ(x − x2 ). Plugging these into (43) gives
Z
Z
δ(x − x1 )G(x; x2 )dx −
G(x; x1 )δ(x − x2 )dx = 0
(45)
Ω
Ω
Using the basic property of the delta function, this simplifies to
G(x1 ; x2 ) − G(x2 ; x1 ) = 0
(46)
which is what we wanted to show.
A related fact has to do with interchanging partial derivatives of G.
Take the one-dimensional case, and observe
∂x G(x, x0 ) = lim (G(x + h, x0 ) − G(x, x0 ))/h = lim (G(x0 , x + h) − G(x0 , x))/h
h→0
h→0
=∂x0 G(x0 , x).
(47)
In other words, reciprocity implies that we can interchange the partial derivatives with respect to x and x0 provided we also interchange the arguments
x and x0 .
14
4.3
Green’s formula representation for inhomogeneous boundary conditions
Let’s now use (44) to solve a problem of the form
Lu(x) = f (x)
+ inhomogeneous boundary conditions.
(48)
We do this by setting v(x) = G(x; x0 ) in the Green’s formula (44), giving
Z
Z
(49)
(LG)u dx − (Lu)G dx = boundary terms.
Ω
Ω
Because LG = δ(x − x0 ) and Lu = f , we have
Z
Z
δ(x − x0 )u(x)dx −
f (x)G(x; x0 )dx = boundary terms
Ω
(50)
Ω
We can collapse the integral involving the δ function, leading to
Z
u(x0 ) =
G(x; x0 )f (x)dx + boundary terms
(51)
Ω
Provided we can evaluate everything on the right hand side, this gives a
formula for the solution which has two components: one which accounts
for the inhomogeneous term in the equation, and another (boundary integral) which accounts for the inhomogeneous boundary conditions.
4.4
Green’s formula representation for the Laplace equation
The discussion above was meant to sketch the structure of the Green’s function solution for a general linear equation with general inhomogeneous
boundary conditions. Now let’s see exactly what this looks like for the
particular case L = ∆.
One dimension. In this case, the Green’s formula (44) is nothing more than
integration by parts twice (i.e. Green’s second identity in 1-D), which for
L = d2 /dx2 and Ω = [0, L] reads
Z
L
uv 00 − vu00 dx = [uv 0 − vu0 ]L
0.
(52)
0
If we want to solve
uxx = f,
u(0) = A,
15
u(L) = B,
(53)
we use the Green’s function (24) appropriate for homogeneous boundary
conditions. Plugging v(x) = G(x; x0 ) into (52), we get
Z
L
u(x)Gxx (x; x0 ) − G(x; x0 )u00 (x)dx = [u(x)Gx (x; x0 ) − G(x; x0 )u0 (x)]x=L
x=0 .
0
(54)
Using uxx = f , Gxx (x; x0 ) = δ(x − x0 ) and noting that G = 0 on the boundaries, this can be written
Z L
u(x0 ) =
G(x; x0 )f (x)dx + [u(x)Gx (x; x0 )]x=L
(55)
x=0 .
0
The first term on the right looks like the formula we had for homogeneous
boundary conditions, with an important exception: x and x0 are in the
wrong places. To make this look more like the original formula (4), we
can interchange the notation for x and x0 ,
L
Z
u(x) =
0
G(x0 ; x)f (x0 )dx0 + [u(x0 )Gx (x0 ; x)]xx00 =L
=0 ,
(56)
(notice the partial derivative on G does not change, since it is just the
derivative with respect to the first input variable regardless of what we
call it!) We can also use the reciprocity identities so that the first term looks
exactly like (4), giving
Z
L
G(x; x0 )f (x0 )dx0 + BGx0 (x; L) − AGx0 (x; 0).
u(x) =
(57)
0
Note this step was for aesthetic purposes only; formula (55) is sufficient to
solve the problem, but (57) is easier to interpret. For example, the influence
of the right hand boundary condition at a point x is Gx0 (x; L).
Many dimensions. In this case, the Green’s formula (44) is nothing more
than Green’s identity
Z
Z
u∆v − v∆u dx =
u∇v · n
ˆ − v∇u · n
ˆ dx.
(58)
Ω
∂Ω
Suppose we wish to solve the problem with the inhomogeneous boundary
condition
∆u = f in Ω, u = h on ∂Ω,
Let G be the Green’s function that solves ∆G = δ(x − x0 ) with homogeneous, Dirichlet boundary conditions. Then substituting v(x) = G(x, x0 )
16
in (58), we have (being careful to keep x as the variable of integration)
Z
Z
u(x)δ(x−x0 )−G(x; x0 )f (x) dx =
u(x)∇x G(x; x0 )·ˆ
n(x)−G(x; x0 )∇u(x)·ˆ
n(x) dx.
Ω
∂Ω
(59)
The notation n
ˆ (x) reminds the reader the unit normal vector is located at
x, not x0 . Using the definition of the δ-function and the fact that G is zero
on the domain boundary, this simplifies to
Z
Z
h(x)∇x G(x; x0 ) · n
ˆ (x) dx.
(60)
G(x; x0 )f (x)dx +
u(x0 ) =
∂Ω
Ω
Again, we can make this look like (4) by interchanging x and x0 and using
reciprocity:
Z
Z
ˆ (x0 ) dx0 .
(61)
h(x0 )∇x0 G(x; x0 ) · n
G(x; x0 )f (x0 )dx0 +
u(x) =
∂Ω
Ω
The first term is just the solution we expect for homogeneous boundary
conditions. The second term is more surprising: it is the derivative of G that
goes into the formula to account for the inhomogeneous Dirichlet boundary condition.
Example: the Poisson integral formula revisited. In the case that Ω is a
disk of radius a, we have found the Green’s function which is zero on the
boundary in equation (41). The boundary value problem
∆u = 0,
u(a, θ) = h(θ)
has a solution given by (61), where f = 0. The normal derivative to the
boundary of G is just the radial derivative
∇x0 G(x; x0 )·ˆ
n(x0 ) = Gr0 (r, θ; r0 , θ0 )
1
2r0 − 2r cos(θ − θ0 )
2r0 r2 − 2ra2 cos(θ − θ0 )
=
−
,
4π r2 + r02 − 2rr0 cos(θ − θ0 ) r2 r02 + a4 − 2rr0 a2 cos(θ − θ0 )
which at r0 = a is
a
2π
1 − (r/a)2
r2 + a2 − 2ar cos(θ − θ0 )
.
The integral on the boundary can be be written explicitly in terms of the
parameter θ0 , and formula (61) becomes
Z 2π
1
(a2 − r2 )h(θ0 )
u(r, θ) =
dθ0 .
2π 0 a2 + r2 − 2ar cos(θ − θ0 )
which is the same formula as we obtained from separation of variables.
17
4.5
Neumann boundary conditions
Suppose we want to solve Laplace’s equation in the upper half space {(x, y, z)|z >
0}, with both a far-field boundary condition and a Neumann condition on
the xy-plane:
∆u = 0,
lim u(x, y, z) = 0,
z→∞
uz (x, y, 0) = h(x, y).
(62)
Notice that the Green’s formula (58) has boundary terms than involve both
Dirichlet and Neumann data. Of course, we only know the derivative of u
on the boundary, so we need to make sure that boundary terms involving
Dirichlet data will vanish. To make this happen, the Green’s function substituted in for v must have ∇G · n
ˆ = 0 on the boundary. In other words,
to use Green’s formulas and functions in conjunction we must respect the
“boundary condition principle”:
The Green’s function must have the same type of boundary conditions
as the problem to be solved, and they must be homogeneous.
For (62), we need a Green’s function which has Gz (x, y, 0; x0 , y0 , z0 ) = 0,
and vanishes at infinity. The method of images tells us this can be done
using the even reflection across the xy-plane. In other words, we want
G(x, y, z; x0 , y0 , z0 ) = G3 (x, y, z; x0 , y0 , z0 ) + G3 (x, y, z; x0 , y0 , −z0 )
1
=
4π
1
1
p
+p
(x − x0 )2 + (y − y0 )2 + (z − z0 )2
(x − x0 )2 + (y − y0 )2 + (z + z0 )2
It is easy to check Gz = 0 when z = 0. Now substituting v(x) = G(x, x0 )
into (58), we obtain after collapsing the δ-function integral
Z
u(x0 ) = −
G(x; x0 )∇u(x) · n
ˆ (x) dx,
(63)
∂Ω
since ∆u = 0 and derivatives of G vanish on the boundary. We should
interpret ∂Ω as both the xy-plane and the effective boundary at infinity, but
on the latter part u and G vanish. Notice that since n
ˆ is directed outward
relative to Ω, ∇u(x)·ˆ
n(x) = −uz (x, y, 0), so that in coordinates (63) becomes
(after interchanging notation between x and x0 and using reciprocity)
Z ∞Z ∞
1
h(x0 , y0 )
p
u(x, y, z) =
dx0 dy0 .
(64)
2π −∞ −∞ (x − x0 )2 + (y − y0 )2 + z 2
This is a type of Poisson formula since it solves Laplace’s equation by a
boundary integral.
18
!
.