## Search in the document preview

**ME201/MTH281/ME400/CHE400
Legendre Polynomials
**

**1. Introduction
**

This notebook has three objectives: (1) to summarize some useful information about Legendre polynomials, (2) to show
how to use *Mathematica* in calculations with Legendre polynomials, and (3) to present some examples of the use of Legendre
polynomials in the solution of Laplace's equation in spherical coordinates. In our course, the Legendre polynomials arose from
separation of variables for the Laplace equation in spherical coordinates, so we begin there. The basic spherical coordinate
system is shown below. The location of a point P is specified by the distance *r* of the point from the origin, the angle f between
the position vector and the *z*-axis, and the angle q from the *x*-axis to the projection of the position vector onto the *xy* plane.

The Laplace equation for a function F(*r*, f, q) is given by

(1)“2F = 1

*r*2
¶∂

¶∂r
*r*2
¶∂F

¶∂r +

1

*r*2 sinf

¶∂

¶∂f sinf

¶∂F

¶∂f +

1

*r*2 sin2 f

¶∂2F

¶∂q2 = 0 .

In this notebook, we will consider only axisymmetric solutions of (1) -- that is, solutions which depend on *r* and f but not on q.
Then equation (1) reduces to

(2) 1

*r*2
¶∂

¶∂r
*r*2
¶∂F

¶∂r +

1

*r*2 sinf

¶∂

¶∂f sinf

¶∂F

¶∂f = 0 .

As we showed in class by a rather lengthy analysis, equation (2) has separated solutions of the form

(3)*rn Pn *HcosfL and *r*-H*n*+1L *Pn *HcosfL ,
where *n *is a non-negative integer and *Pn *is the *n*th Legendre polynomial. These solutions can be used to solve axisymmetric
problems inside a sphere, exterior to a sphere, or in the region between concentric spheres. We include examples of each type
later in this notebook. Now we look in more detail at Legendre's equation and the Legendre polynomials.

**2. Legendre Polynomials
**

**ü 2.1 Differential Equation
**

The first result in the search for separated solutions of equation (2), which ultimately leads to the formulas (3), is the pair
of differential equations (4) for the *r*-dependent part *F*(*r*), and the f-dependent part *P*(f) of the separated solutions:

(4)*r*2
d2 *F
*

d*r*2
+ 2 *r
*

d*F
*

d*r
*- l *F *= 0 , and

d

df sinf

d*P
*

df
+ lsinf *P *= 0 ,

where l is the separation constant. The *r*-equation is equidimensional and thus has solutions, easily found, which are powers of
*r*. The f-equation is Legendre's equation. We begin by transforming it to a somewhat simpler form by a change of independent
variable, namely

(5)h = cosf .

Then equation (4) becomes

(6) d

dh
BI1 - h2M d*P
*

dh
F + l *P *= 0 , -1 < h < 1 .

Equation (6) has regular singular points at the endpoints h = ± 1, and we require the solution *P* to be well-behaved at those
points. This is our first example of a singular Sturm-Liouville system. Those special values of l for which there are such well-
behaved solutions are the eigenvalues of the problem.

As we showed in class, we may find solutions of (6) in the form of power series about h = 0: *P*(h) = ⁄*n*=0¶ *an *h*n*. By
substitution of this into equation (6), we find the recurrence relation for the coefficients *an *:

(7)*an*+2 =
*n *H*n *+ 1L - l
H*n *+ 1L H*n *+ 2L *an *.

Because of the index increment of 2 in (7), the solutions fall naturally into even and odd functions of h, with the coefficient
sequences of the two not mixing. For very large *n*, we may approximate the relation (7) by ignoring l compared with *n*(*n* + 1).
The result is (*n* + 2)*an*+2> *nan *-- that is, *an*> constant/*n*. From that result it is possible to show that any such solution is logarithmi-
cally singular at one or both endpoints. This result is suggested by (but not proved by) the two series

(8)ln H1 - hL = -‚ n=1

¶ h*n
*

*n
*, and ln H1 + hL = ‚

n=1

¶ H-1L*n*+1 h*n
n
*

.

The only way to avoid such singularities in the solution is for the series to terminate. We see immediately from the recurrence
relation (7) that termination occurs if and only if l = *k(k* + 1) for *k* a non-negative integer. The terminating solution in that case
is a polynomial of degree *k*. If *k* is even, the polynomial has only even powers and is then an even function of h. If *k* is odd, only
odd powers appear and the function is odd. Such solutions are called Legendre polynomials. The *k*th one is denoted by *Pk*(h),
with the convention that the arbitrary multiplicative constant is fixed by the condition

(9)*Pk *H1L = 1 .
Any of the polynomials can be constructed directly from the recurrence formula (7) and the normalization (9). As an alternative,
there is the well-known formula of Rodrigues, which gives an explicit expression for the *n*th polynomial. Although it is not
usually used to compute the polynomials, it is still of interest:

**2 *** legendre.nb*

Any of the polynomials can be constructed directly from the recurrence formula (7) and the normalization (9). As an alternative,
there is the well-known formula of Rodrigues, which gives an explicit expression for the *n*th polynomial. Although it is not
usually used to compute the polynomials, it is still of interest:

(10)*Pk *HhL = 1
2*k k *!

d*k
*

dh*k
*Ih2 - 1M*k *.

The Legendre polynomials are built into *Mathematica*. *Mathematica*'s notation is LegendreP[k,h] for *Pk*(h). We now use
*Mathematica* to obtain the formulas for the first 11 of these polynomials. We put them in a table.

**TableForm@Table@8i, i * Hi + 1L, LegendreP@i, hD<, 8i, 0, 10<D,
TableHeadings Ø 8None, 8"k", "lk", "PkHhL"<<D
**

k lk PkHhL 0 0 1 1 2 h

2 6 1 2 I-1 + 3 h2M

3 12 1 2 I-3 h + 5 h3M

4 20 1 8 I3 - 30 h2 + 35 h4M

5 30 1 8 I15 h - 70 h3 + 63 h5M

6 42 1 16

I-5 + 105 h2 - 315 h4 + 231 h6M 7 56 1

16 I-35 h + 315 h3 - 693 h5 + 429 h7M

8 72 1 128

I35 - 1260 h2 + 6930 h4 - 12 012 h6 + 6435 h8M 9 90 1

128 I315 h - 4620 h3 + 18 018 h5 - 25 740 h7 + 12 155 h9M

10 110 1 256

I-63 + 3465 h2 - 30 030 h4 + 90 090 h6 - 109 395 h8 + 46 189 h10M By inverting these we get a useful set of formulas expressing the powers of h in terms of the Legendre polynomials. The

first five are given below.

(11)

1 = *P*0HhL
h = *P*1HhL
h2 =

1

3
H*P*0HhL + 2 *P*2HhLL

h3 = 1

5
H3 *P*1HhL + 2 *P*3HhLL

h4 = 1

35
H7 *P*0HhL + 20 *P*2HhL + 8 *P*4HhLL

These formulas are useful in obtaining Legendre expansions of polynomials. A second equally useful set of formulas may be
derived from these and from the trig addition and double angle formulas. In this second set of formulas, we express cos(*n*f) in
terms of *Pk*(cosf), with *n* an integer and *k* § *n.* The first five of those formulas are given below.

*legendre.nb ***3**

(12)

cos H0 ÿfL = 1 = *P*0HcosfL
cos HfL = *P*1HcosfL
cos H2 fL = 1

3
H-*P*0HcosfL + 4 *P*2HcosfLL

cos H3 fL = 1
5
H-3 *P*1HcosfL + 8 *P*3HcosfLL

cos H4 fL = 1 105

H-7 *P*0HcosfL - 80 *P*2HcosfL + 192 *P*4HcosfLL

These last two sets of formulas are especially useful in solving boundary value problems in spherical coordinates when the boundary conditions have special, simple forms.

**ü 2.2 Some Useful Formulas and Graphs
**

To get some idea of what these polynomials look like, we construct graphs of the first 7. According to the general result about the zeros of solutions of Sturm-Liouville systems, the kth polynomial should have exactly k zeros in the interval (-1,1). We first define a function legraph[k] that produces a graph of the kth polynomial, and then we use a Do loop to construct the first 7 graphs. The size of the graphs throughout this notebook is controlled by the SetOptions command below. The value 250 is appropriate for a printed version. A value of 350 might be better for computer display.

**SetOptions@Plot, ImageSize Ø 250D;
legraph@k_D := Plot@LegendreP@k, hD, 8h, -1, 1<, AxesLabel Ø 8"h", "PkHhL"<,
**

**PlotRange Ø 88-1.1, 1.1<, 8-1.1, 1.1<<, PlotLabel Ø Row@8"k =", PaddedForm@k, 3D<DD
**

**4 *** legendre.nb*

**GraphicsGrid@Table@8legraph@iD, legraph@i + 1D<, 8i, 0, 6, 2<DD
**

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

*Pk*HhLk = 0

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

*Pk*HhLk = 1

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

*Pk*HhLk = 2

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

*Pk*HhLk = 3

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

*Pk*HhLk = 4

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

*Pk*HhLk = 5

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

*Pk*HhLk = 6

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

*Pk*HhLk = 7

We see the expected alternation between even and odd functions, and the expected number of zeros in each case.

There are a large number of formulas involving Legendre polynomials. We consider here only a few of the most useful. The following is a recurrence relation for three consecutive Legendre polynomials:

(13)H*n *+ 1L *Pn*+1 HhL - H2 *n *+ 1L h*Pn *HhL + *nP*n-1 HhL = 0 .
Although *Pn*-1 is not defined when *n* = 0, you can easily check that the formula (13) remains true in that case if we define *P*-1 to
be anything bounded. An interesting application of (13) is to compute *Pn*HhL, starting with the known functions *P*0HhL = 1 and
*P*1HhL = h. Then we get

*legendre.nb ***5**

Although *Pn*-1 is not defined when *n* = 0, you can easily check that the formula (13) remains true in that case if we define *P*-1 to
be anything bounded. An interesting application of (13) is to compute *Pn*HhL, starting with the known functions *P*0HhL = 1 and
*P*1HhL = h. Then we get

(14)

*Pn*+1 HhL = H2 *n *+ 1L h *Pn *HhL - *n Pn*-1 HhL
*n *+ 1

so *P*2 HhL = 3 h *P*1 HhL - *P*0 HhL
2

= 3 h2 - 1

2 ,

and *P*3 HhL = 5 h *P*2 HhL - 2 *P*1 HhL
3

= 5

3 h

3 h2 - 1

2 -

2

3 h =

5

2 h3 -

3

2 h .

This can be continued to any order.

A second useful recurrence relation for three consecutive functions contains derivatives:

(15)*P*£*n*+1 HhL - *P*£*n*-1 HhL = H2 *n *+ 1L *Pn *HhL .
This formula is valid for *n* r 1.

Now we look at values of the polynomials for special values of h. Because of the way that the polynomials are defined, we have

(16)*Pn *H1L = 1 .
*Pn* is an even function for *n* even and an odd function for *n* odd, so

(17)*Pn *H-1L = H-1Ln .
Now consider h = 0. When *n* is odd, *Pn* is odd, so that *Pn*(0) = 0. When *n* is even, we may apply (13) repeatedly, using h = 0,
and starting with *n* = 1, to get

(18)*P*2 *n *H0L = H-1L
*n *H2 *n *L !

22 *n *H*n *!L2 .
We use these results to calculate a few integrals involving Legendre polynomials. These integrals will be useful in

constructing expansions in the next section. We first use (15) and (17) to construct an indefinite integral of *Pn*, starting from a
lower limit of -1. We get

(19)‡ -1

h

*Pn *HhèL „hè = *Pn*+1 HhL - *Pn*-1 HhL
2 *n *+ 1

.

This is valid for *n* r 1. It follows immediately from (16) and (19) that

(20)‡ -1

1
*Pn *HhL „h = 0 ,

for *n* ¥ 1. This is a special case of the orthogonality discussed in the next section. (Because *P*0= 1, *Pn* = *PnP*0, and (17) then says
that *P*0 and *Pn* are orthogonal on[-1,1].) Another special case of (19) is h = 0, which gives

(21)‡ -1

0
*Pn *HhL „h = *Pn*+1 H0L - *Pn*-1 H0L

2 *n *+ 1
.

If *n* is even the right-hand side is zero. If *n* is odd, the values at 0 are given by (18). By subtracting (21) from (20) we get

**6 *** legendre.nb*

(22)‡ 0

1
*Pn *HhL „h = *Pn*-1 H0L - *Pn*+1 H0L

2 *n *+ 1
,

which again is zero for *n* even, and can be evaluated for *n* odd by using (18).

We can also evaluate integrals of the form Ÿ hP*n*HhL „h by using the recurrence relations (13) and (15). The result of the
somewhat tedious calculation is the formula below, valid for *n* r 2.

(23)‡ -1

h

hè *Pn *HhèL „hè = H*n *+1L H2 *n *-1L *Pn*+2 HhL + H2 *n *+ 1L *Pn *HhL - *n *H*2 n *+3L *Pn*-2 HhLH2 *n *-1L H2 *n *+ 1L H2 *n *+3L .
For the special case h = 1, we get from (16) and (23)

(24)‡ -1

1
h*Pn *HhL „h = 0 ,

valid for *n* ¥ 2, again a special case of orthogonality, because *P*1HhL = h. Another special case of (23) is h = 0, which gives
(25)‡

-1

0
h*Pn *HhL „h = H*n *+1L H2 *n *-1L *Pn*+2 H0L + H2 *n *+ 1L *Pn *H0L - *n *H*2 n *+3L *Pn*-2 H0LH2 *n *-1L H2 *n *+ 1L H2 *n *+3L ,

which is zero for *n* odd, and may be evaluted from (18) for *n* even. Finally, by subtracting (25) from (24), we get

(26)‡ 0

1
h*Pn *HhL „h = *n *H*2 n *+3L *Pn*-2 H0L - H2 *n *+ 1L *Pn *H0L - H*n *+1L H2 *n *-1L *Pn*+2 H0LH2 *n *-1L H2 *n *+ 1L H2 *n *+3L .

We will use some of these integrals in the examples of expansions in the next section, and in the Laplace equation examples in the last section.

**ü 2.3 Expansions in Legendre Polynomials
**

As we showed in class from the differential equation (6), the Legendre polynomials are orthogonal on the interval [-1,1]:

(27)‡ -1

1
*Pm *HhL *Pn *HhL „h = 0 for *m *¹≠ *n *.

It may be shown that the normalization integral is given by

(28)‡ -1

1
*Pn *HhL *Pn *HhL „h = 2

2 *n *+ 1
.

The polynomials form a complete set on the interval [-1,1], and any piecewise smooth function may be expanded in a series of the polynomials. The series will converge at each point to the usual mean of the right and left limits. The coefficients are easily calculated using the orthogonality properties (24) and the normalization integral (25), and we have

(29)*f *HhL = ‚
n=0

¶

*Cn Pn *HhL , where *Cn *= 2 *n *+ 1
2

‡ -1

1
*f *HhL *Pn *HhL „h .

As an example, we expand the step function given below in such a series.

(30)*f *HhL = -1 for - 1 b h b 0, and *f *HhL = 1 for 0 < h § 1.
Because *f* is odd, the even coefficients will all vanish. The odd coefficients are given by

(31)*Cn *= H 2 *n *+ 1L ‡
0

1
*Pn *HhL „h .

By using equation (22) we get

*legendre.nb ***7**

By using equation (22) we get

(32)*Cn *= @*Pn*-1 H0L - *Pn*+1 H0LD ,
where the values at h = 0 needed in (32) are given by equation (18). We now use *Mathematica* to calculate the coefficients up to
*n *= 51, starting with *n* = 1.

**Module@8i<, Co = 8<; Do@Co = Append@Co,HLegendreP@i - 1, 0D - LegendreP@i + 1, 0DLD, 8i, 1, 51<DD
**The kth coefficient is then given by Co[[k]]. We sample a few values:

**Co@@1DD
**3

2

**Co@@2DD
**0

**Co@@9DD
**133

256

All the even coefficients are zero.

Now we define the kth partial sum, and then a graph of the kth partial sum, along with the original function. Because the
coefficient of *P*0 is zero, we may start the sum over* i* with the *i* = 1 term. The given function *f*(h) is shown in blue and the partial
sums in red.

**legsum@h_, k_D := Module@8i<, Sum@N@Co@@iDDD * LegendreP@i, hD, 8i, 1, k<DD
f@h_D := If@Hh < 0L, H-1L, H1LD
legraph@k_D := Plot@8f@hD, legsum@h, kD<, 8h, -1, 1<, PlotRange Ø 8-1.5, 1.5<,
**

**PlotStyle Ø 88RGBColor@0, 0, 1D, Thickness@0.004D<, 8RGBColor@1, 0, 0D, Thickness@0.004D<<,
AxesLabel Ø 8"h", "fHhL"<, PlotLabel Ø Row@8"k =", PaddedForm@k, 2D<DD
**

We construct a sequence of partial sums which may be animated to see the convergence. We go up to *P*51. The even terms are
zero, so we increment by 2 in the sequence of partial sums. It would be much more efficient computationally to save each partial
sum, and then use it to compute the next partial sum in the sequence. The present inefficient method is however much easier to
program. The output from the calculation is sent to a Manipulate panel.

**8 *** legendre.nb*

**DynamicModule@8mangraph, i<, Do@mangraph@iD = legraph@iD, 8i, 1, 51, 2<D;
Manipulate@mangraph@iD, 8i, 1, 51, 2<DD
**

i

-1.0 -0.5 0.5 1.0 h

-1.5

-1.0

-0.5

0.5

1.0

1.5 fHhLk = 1

You can use the slider to move through the graphs, and also to display a movie of the graphs. When we do this, we see the familiar Gibbs overshoot at the discontinuity. We also see a struggle to converge at the endpoints.

For visualization in the printed notebook, we print out every 10th graph in the sequence.

*legendre.nb ***9**

**GraphicsGrid@Table@8legraph@iD, legraph@i + 10D<, 8i, 1, 51, 20<DD
**

-1.0 -0.5 0.5 1.0 h

-1.5

-1.0

-0.5

0.5

1.0

1.5 fHhLk = 1

-1.0 -0.5 0.5 1.0 h

-1.5

-1.0

-0.5

0.5

1.0

1.5 fHhLk = 11

-1.0 -0.5 0.5 1.0 h

-1.5

-1.0

-0.5

0.5

1.0

1.5 fHhLk = 21

-1.0 -0.5 0.5 1.0 h

-1.5

-1.0

-0.5

0.5

1.0

1.5 fHhLk = 31

-1.0 -0.5 0.5 1.0 h

-1.5

-1.0

-0.5

0.5

1.0

1.5 fHhLk = 41

-1.0 -0.5 0.5 1.0 h

-1.5

-1.0

-0.5

0.5

1.0

1.5 fHhLk = 51

As a second example, we consider a function which is continuous on the interval [-1,1], but which has a discontinuity in slope at 0. The function is defined to be 0 in the left-half interval [-1,0] and h in the right half interval:

**y@h_D := If@Hh < 0L, H0L, HhLD
**The expansion coefficients are given by

(33)*Dn *=
H 2 *n *+ 1L

2 ‡ 0

1
h*Pn *HhL „h .

By equation (26) we get

(34)*Dn *=
*n *H*2 n *+3L *Pn*-2 H0L - H2 *n *+ 1L *Pn *H0L - H*n *+1L H2 *n *-1L *Pn*+2 H0L

2 H2 *n *-1L H2 *n *+3L .
This is valid for *n* r 2. The first two coefficients are

(35)

**10 *** legendre.nb*

(35)

*D*0 =
1

2 ‡ 0

1
h*P*0 HhL „h = 14 ,

and *D*1 =
3

2 ‡ 0

1
h*P*1 HhL „h = 12 .

We now construct an array containing the first 20 coefficients. For convenience, we give a name to the right-hand side of (34).

**coeff@n_D :=HHn H2 n + 3L LegendreP@n - 2, 0D - H2 n + 1L LegendreP@n, 0DL - Hn + 1L H2 n -1L LegendreP@n + 2, 0DL êH2 H2 n -1L H2 n +3LL
**We construct the first 20 coefficients, starting at *n* = 1, and storing them in the array Co2.

**Module@8i<, Co2 = 80.5<; Do@Co2 = Append@Co2,
N@coeff@iDDD, 8i, 2, 20<DD
**

The module below defines the kth partial sum of the Legendre expansion and assigns it to legsum2. The 0.25 appearing in the
module is the *n* = 0 term added to the sum.

**legsum2@h_, k_D := Module@8i<, 0.25 + Sum@N@Co2@@iDDD * LegendreP@i, hD, 8i, 1, k<DD
**We define a function graph[k] which graphs the kth partial sum in red and the original function in blue. Then we use a Do loop
to construct a sequence of partial sums. All of the graphs of the sequence are collected in a Manipulate panel.

**graph@k_D := Plot@8y@hD, legsum2@h, kD<, 8h, -1, 1<, PlotRange Ø 8-1.1, 1.1<,
PlotStyle Ø 88RGBColor@0, 0, 1D, Thickness@0.004D<, 8RGBColor@1, 0, 0D, Thickness@0.004D<<,
AxesLabel Ø 8"h", "yHhL"<, PlotLabel Ø Row@8"k =", PaddedForm@k, 2D<DD
**

**DynamicModule@8mangraph, i<, Do@mangraph@iD = graph@iD, 8i, 0, 10, 1<D;
Manipulate@mangraph@iD, 8i, 0, 10, 1<DD
**

i

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

yHhLk = 0

We see that in this case, the convergence is very rapid. For convergence to graphical accuracy, 10 terms are sufficient. There is far less flailing around at the endpoints than in the previous case. For visualization in the printed version of this note- book, we print out every other graph.

*legendre.nb ***11**

**GraphicsGrid@Table@8graph@iD, graph@i + 2D<, 8i, 0, 8, 4<DD
**

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

yHhLk = 0

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

yHhLk = 2

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

yHhLk = 4

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

yHhLk = 6

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

yHhLk = 8

-1.0 -0.5 0.5 1.0 h

-1.0

-0.5

0.5

1.0

yHhLk = 10

In both of our examples, we have used the full machinary of the orthogonality and integration to calculated the coeffi-
cients. In cases where *f*(h) is a polynomial, we may take advantage of the formulas (11) to obtain the desired expansions. For
example, if

(36)*f *HhL = h4,
then we get directly from the last of equations (11)

(37)*f *HhL = 1
35

H7 *P*0 HhL + 20 *P*2 HhL + 8 *P*4 HhLL.
Throughout this section, we have used the modified independent variable h, defined by equation (5). The expansion

theorem in terms of this variable is

(38)*f *HhL = ‚
n=0

¶

*Cn Pn *HhL , where *Cn *= H 2 *n *+ 1L
2

‡ -1

1
*f *HhL *Pn *HhL „h .

If we recast this in terms of the original variable f, where h = cosf and *f*(h) = *f*(cosf) = *g*(f), we get

**12 *** legendre.nb*

If we recast this in terms of the original variable f, where h = cosf and *f*(h) = *f*(cosf) = *g*(f), we get

(39)*g *HfL = ‚
n=0

¶

*Cn Pn *HcosfL , where *Cn *= H 2 *n *+ 1L
2

‡ 0

p

*g *HfL *Pn *HcosfL sinf „f .
As the notation suggests, the coefficients *Cn* are the same in the expansions (38) and (39).

We will now use these expansions in solving the Laplace equation.

**3. Examples of Solutions of Laplace's Equation
**

**ü 3.1 Interior of a Sphere
**

For our first example, we find the electrostatic potential F inside a sphere of radius *a*, when the potential on the boundary
is a particularly simple given function. The problem statement is given below.

(40)

1

*r*2
¶∂

¶∂r
*r*2
¶∂F

¶∂r +

1

*r*2 sinf

¶∂

¶∂f sinf

¶∂F

¶∂f
= 0 , *r *< *a*, and 0 § f § p ,

with F H*a*, fL = *g *HfL = F0 H1 + 2 cos HfL + 3 cos H2 fL L,
where F0is a constant. This is a situation in which the formulas (13) are particularly useful. With those formulas, we may
express each of the three angular harmonics in (40) in terms of the Legendre polynomials with argument cosf. The result is

(41)*g *HfL = F0 H2 *P*1 HcosfL + 4 *P*2 HcosfLL.
The separated solutions are given by equation (3). The solutions with negative powers of *r* are singular at the origin and are thus
not suitable for the potential inside the sphere. We superpose all of the solutions with non-negative powers to get F:

(42)F H*r*, fL = ‚
n=0

¶

*An rn Pn *HcosfL .
We now imposed the boundary condition:

(43)F0 H2 *P*1 HcosfL + 4 *P*2 HcosfLL = ‚
n=0

¶

*A*n *a*n *P*n HcosfL,
which gives

(44)*A*1 = 2 F0 ê*a*, *A*2 = 4 F0 ë*a*2, and all other *A*n = 0.
The solution is then

(45)F H*r*, fL = 2 F0AH*r *ê*a*L *P*1 HcosfL + 2 H*r *ê*a*L2 *P*2 HcosfL.

*legendre.nb ***13**

For our second example, we find the electrostatic potential F inside a sphere of radius *a*, when the potential on the
boundary is the step function *g*(f) = *V*0 for 0 bf b p/2, and -*V*0 for p/2 < f b p. The problem statement is

(46)

1

*r*2
¶∂

¶∂r
*r*2
¶∂F

¶∂r +

1

*r*2 sinf

¶∂

¶∂f sinf

¶∂F

¶∂f
= 0 , *r *< *a*, and 0 § f § p ,

with F H*a*, fL = *V*0 for 0 § f § p
2

,

and F H*a*, fL = -*V*0 for p
2

< f § p .

The separated solutions are given by equation (3). The solutions with negative powers of *r* are singular at the origin and are thus
not suitable for the potential inside the sphere. We superpose all of the solutions with non-negative powers to get F:

(47)F H*r*, fL = ‚
n=0

¶

*An rn Pn *HcosfL .
Imposing the boundary condition, we get

(48)F H*a*, fL = ‚
n=0

¶

*An an Pn *HcosfL .
Apart from the factor of the voltage *V*0, the boundary function in this case is the function expanded in the first example of section
2.3, and we found there the coefficents *Cn* for the expansion. (The numerical values of the coefficients were stored in the array
Co.) Thus

(49)*An an *= *V*0 *Cn ,
*

where *Cn* is given by equation (32), and is equal to Co[[n]]. Then the kth partial sum of the series for F is given by

**Fsum@r_, f_, k_D := V0 SumAHr ê aLi Co@@iDD LegendreP@i, Cos@fDD, 8i, 1, k<E
**We look at the first few partial sums.

**Fsum@r, f, 1D
**3 r V0 Cos@fD

2 a

**Fsum@r, f, 3D
**V0

3 r Cos@fD 2 a

- 7 r3 I-3 Cos@fD + 5 Cos@fD3M

16 a3

**Fsum@r, f, 5D
**V0

3 r Cos@fD 2 a

- 7 r3 I-3 Cos@fD + 5 Cos@fD3M

16 a3 + 11 r5 I15 Cos@fD - 70 Cos@fD3 + 63 Cos@fD5M

128 a5

In order to evaluate our solution in detail, we choose some specific values for the sphere radius *a*, and the boundary
voltage *V*0. We take

**a = 2.0 H** m **L; V0 = 5.0 H** volts **L;
**First we check our boundary condition. We plot the voltage on the surface *r* = *a* of the sphere. We use terms up to *n *= 51 in the
series, using the coefficients that we calculated earlier. We ask for 200 sample points rather than accepting the default of 25.
The plot takes a very long time because of all the evaluations of both the *Pn*'s and the cosines.

**14 *** legendre.nb*

First we check our boundary condition. We plot the voltage on the surface *r* = *a* of the sphere. We use terms up to *n *= 51 in the
series, using the coefficients that we calculated earlier. We ask for 200 sample points rather than accepting the default of 25.
The plot takes a very long time because of all the evaluations of both the *Pn*'s and the cosines.

**PlotBFsum@a, f, 51D, 8f, 0, p<, PlotRange Ø 8-6, 6<, AxesLabel Ø 8"f", "FHa,fL"<,
PlotPoints Ø 200, Ticks Ø ::0, p
**

**4
,
**

**p
**

**2
,
3 p
**

**4
, p>, Automatic>F
**

p

4

p

2

3 p

4 p

f

-6

-4

-2

2

4

6 FHa,fL

We see that the trend is correct, but that many more terms would be needed for an accurate representation of the boundary
condition. Fortunately, when we evaluate the solution away from the boundary, we have the factor (*r*/*a*L*i* in the *i*th term and this
greatly accelerates the convergence of the series. To show this, we look at the solution on an interior sphere, namely *r* = 1. We
plot a short sequence of partial sums of F versus f for *r* = 1. We use 5, 11, and 17 terms in the partial sums.

**DoBPrintBPlotBFsum@1.0, f, 5 + 6 * iD, 8f, 0, p<, PlotRange Ø 8-4, 4<,
AxesLabel Ø 8"f", "FH1,fL"<, PlotPoints Ø 100, Ticks Ø ::0, p
**

**4
,
**

**p
**

**2
,
3 p
**

**4
, p>, Automatic>,
**

**PlotLabel Ø Row@8"k =", PaddedForm@5 + 6 * i, 3D<DFF, 8i, 0, 2<F
**

*legendre.nb ***15**

p

4

p

2

3 p

4 p

f

-4

-2

2

4 FH1,fL k = 5

p

4

p

2

3 p

4 p

f

-4

-2

2

4 FH1,fL k = 11

p

4

p

2

3 p

4 p

f

-4

-2

2

4 FH1,fL k = 17

The three curves are graphically identical, as you can see by animating the sequence. Thus at this value of *r*, the solution is well-
represented by three nonzero terms.

**ü 3.2 Exterior of a Sphere
**

We consider in this section the solution of Laplace's equation exterior to a sphere of radius *a*, with the boundary condi-
tion on the sphere being the same as the one in the preceeding problem. The full problem statement is given below.

(50)

1

*r*2
¶∂

¶∂r
*r*2
¶∂F

¶∂r +

1

*r*2 sinf

¶∂

¶∂f sinf

¶∂F

¶∂f
= 0 , *r *> *a*, and 0 § f § p ,

with F H*a*, fL = *V*0 for 0 § f § p
2

,

and F H*a*, fL = -*V*0 for p
2

< f § p ,

and F H*r*, fL ö
*r*z¶

0 .

The separated solutions are given by equation (3). Because of the condition at ¶, only the negative powers of *r* may be used in
this region exterior to a sphere. We superpose all of those solutions to get F:

**16 *** legendre.nb*

The separated solutions are given by equation (3). Because of the condition at ¶, only the negative powers of *r* may be used in
this region exterior to a sphere. We superpose all of those solutions to get F:

(51)F H*r*, fL = ‚
n=0

¶

*Bn r*-H*n*+1L *Pn *HcosfL .
Apart from the factor of the voltage *V*0, the boundary function is the function expanded in the first example of section 2.3, and
the coefficients *Cn* were obtained in that section. Then imposing the boundary condition we get

(52)F H*a*, fL = V0‚
n=0

¶

*Cn Pn *HcosfL = ‚
n=0

¶

*Bn a*-H*n*+1L *Pn *HcosfL,
(53)so *Bn a*-H*n*+1L = *V*0 *Cn ,
*

where *Cn* is given by equation (32). The numerical values of the coefficients are stored in the array Co[[n]] (except for the first
coefficient which is zero in this case). Then the kth partial sum of the series for F is given by

**Fsum2@r_, f_, k_D := V0 SumAHa ê rLi+1 Co@@iDD LegendreP@i, Cos@fDD, 8i, 1, k<E
**We look at the first term in the series after clearing the numerical values assigned earlier to *a* and *V*0.

**Clear@a, V0D;
Fsum2@r, f, 1D
**3 a2 V0 Cos@fD

2 r2

For *r* p *a *this is a reasonable approximation to the entire solution, because the omitted terms are higher inverse powers of *r*.
Thus when we are far from the sphere it looks like a dipole. Higher approximations are obtained by keeping more terms.

**Fsum2@r, f, 3D
**V0

3 a2 Cos@fD 2 r2

- 7 a4 I-3 Cos@fD + 5 Cos@fD3M

16 r4

**Fsum2@r, f, 5D
**V0

3 a2 Cos@fD 2 r2

- 7 a4 I-3 Cos@fD + 5 Cos@fD3M

16 r4 + 11 a6 I15 Cos@fD - 70 Cos@fD3 + 63 Cos@fD5M

128 r6

**ü 3.3 The Region Between Two Concentric Spheres
**

As our final example, we consider the region between two concentric spheres, with radii *a* and *b* , *b* > *a*. We solve the
Laplace equation in the region between the spheres, subject to a boundary condition on each sphere. The problem statement is
given below.

(54)

1

*r*2
¶∂

¶∂r
*r*2
¶∂F

¶∂r +

1

*r*2 sinf

¶∂

¶∂f sinf

¶∂F

¶∂f
= 0 , *a *< *r *< *b*, and 0 § f § p ,

with F H*a*, fL = *g *HfL and F Hb, fL = *h *HfL.
We will specify *g* and *h* explicitly shortly. We leave them general now because that makes it easier to see the structure of the
calculation that way. We begin by expanding both *g* and *h* in Legendre polynomials. That will make our task easier.

*legendre.nb ***17**

(55)*g *HfL = ‚
n=0

¶

*Cn Pn *HcosfL , where *Cn *= H 2 *n *+ 1L
2

‡ 0

p

*g *HfL *Pn *HcosfL sinf „f ,
and

(56)*h *HfL = ‚
n=0

¶

*Dn Pn *HcosfL , where *Dn *= H 2 *n *+ 1L
2

‡ 0

p

*h *HfL *Pn *HcosfL sinf „f .
The coefficients *Cn*and *Dn*are known from the known boundary functions *g* and *h*. The separated solutions are given by equation
(3). The domain of the present problem, *a *<* r *< *b*, does not include either the origin or the point at infinity. Thus there are no
grounds for discarding any of the solutions and we keep them all. The solution for F is then obtained by superposition:

(57)F H*r*, fL = ‚
n=0

¶ I*An rn *+ *Bn r*-H*n*+1LM *Pn *HcosfL .
We now impose the two boundary conditions. We get

(58)F H*a*, fL = ‚
n=0

¶ I*An an *+ *Bn a*-H*n*+1LM *Pn *HcosfL = *g *HfL = ‚
n=0

¶

*Cn Pn *HcosfL ,
and

(59)F H*b*, fL = ‚
n=0

¶ I*An bn *+ *Bn b*-H*n*+1LM *Pn *HcosfL = *h *HfL = ‚
n=0

¶

*Dn Pn *HcosfL .
By equating the coefficients of corresponding terms in the Legendre expansions of (47) and (48), we get

(60)*An an *+ *Bn a*-H*n*+1L = *Cn *, and *An bn *+ *Bn b*-H*n*+1L = *Dn *.
These constitute two linear algebraic equations for each pair *An* and *Bn*. We solve them to get

(61)*An *=
*bn*+1 *Dn *- *an*+1 *Cn
*

*b*2 *n*+1 - *a*2 *n*+1
, and *Bn *=

H*ab*L*n*+1 H*bn Cn *- *an Dn*L
*b*2 *n*+1 - *a*2 *n*+1

.

Now *An* and *Bn* are expressed in terms of known quantities and the solution is complete.

We look at a specific example. We take the functions given below for *g* and *h*.

(62)

*g *HfL = *Va *cos HfL for 0 § f § p
2

,

and *g *HfL = 0 for p
2

< f § p ,

and *h *HfL = 0 for all f.
This means that

(63)*Cn *= *Va *Co2@@nDD *, Dn *= 0 ,
where Co2[[n]] are the coefficients calculated for the second example in section 2.3. We will use the coefficients in a somewhat
different way here. Our new way will be more convenient with respect to the indexing, in that it will allow us to use index zero
for the first term in the series. In section 2.3, we defined a function coeff[n] which returned the value of the nth coefficient
Co2[[n]] in the above expansion. We use that function now to creat the coefficient functions A[n] and B[n] corresponding to
*An *and *Bn*.

**18 *** legendre.nb*

**A@n_D := Va I-an+1 coeff@nD ë Ib2 n + 1 - a2 n + 1MM
B@n_D := Va Ian+1 b2 n + 1 coeff@nD ë Ib2 n + 1 - a2 n + 1MM
**We now define the kth partial sum of the solution F. We call it Fsum3.

**Fsum3@r_, f_, k_D := SumAIA@iD ri + B@iD r-Ii + 1MM LegendreP@i, Cos@fDD, 8i, 0, k<E
**We check this by looking at the first two partial sums.

**Fsum3@r, f, 1D
**-

a Va

4 H-a + bL + a b Va

4 H-a + bL r + a2 b3 Va

2 I-a3 + b3M r2 - a2 r Va

2 I-a3 + b3M Cos@fD
**Fsum3@r, f, 2D
**-

a Va

4 H-a + bL + a b Va

4 H-a + bL r + a2 b3 Va

2 I-a3 + b3M r2 - a2 r Va

2 I-a3 + b3M Cos@fD + 1

2

5 a3 b5 Va

16 I-a5 + b5M r3 - 5 a3 r2 Va

16 I-a5 + b5M I-1 + 3 Cos@fD 2M

Now we will use our solution to calculate the potential on a sphere half way between our two boundary spheres. In
addition, we will check our boundary conditions on *r* = *a* and *r* = *b*. We assign numerical values to the parameters:

**a = 2 H** m **L; b = 4 H** m **L; Va = 10.0 H** volts **L;
**Next we define a function grapher[r,k] which uses the kth partial sum to construct a plot of potential versus f on the sphere of
radius r.

**grapher@r_, k_D := PlotBFsum3@r, f, kD, 8f, 0, p<,
AxesLabel Ø 8"f", "FHr,fL"<, PlotLabel Ø Row@8"r =", PaddedForm@N@rD, 84, 2<D<D,
PlotRange Ø 80, 1.1 * Va<, Ticks Ø ::0, p
**

**4
,
**

**p
**

**2
,
3 p
**

**4
, p>, Automatic>F
**

Now we use 10 terms in the series, and we construct a sequence of 6 plots going in equal *r*-increments from *r* = *a* to *r* = *b*.

*legendre.nb ***19**

**GraphicsGrid@
Table@8grapher@a + i * Hb - aL ê 5, 10D, grapher@a + Hi + 1L * Hb - aL ê 5, 10D<, 8i, 0, 5, 2<DD
**

0 p 4

p

2

3 p

4 p

f

2

4

6

8

10

FHr,fL r = 2.00

0 p 4

p

2

3 p

4 p

f

2

4

6

8

10

FHr,fL r = 2.40

0 p 4

p

2

3 p

4 p

f

2

4

6

8

10

FHr,fL r = 2.80

0 p 4

p

2

3 p

4 p

f

2

4

6

8

10

FHr,fL r = 3.20

0 p 4

p

2

3 p

4 p

f

2

4

6

8

10

FHr,fL r = 3.60

0 p 4

p

2

3 p

4 p

f

2

4

6

8

10

FHr,fL r = 4.00

The first and last graphs verify the boundary conditions that we have imposed on the inner and outer sphere. The remaining graphs show how the solution of the Laplace equation interpolates smoothly between these. We can also use the Manipulate command. We will construct 21 graphs with r varying in equal increments from the inner to the outer boundary.

**20 *** legendre.nb*

**DynamicModule@8mangraph, i, rarg<,
Do@rarg = a + HHb - aL ê 20L * i; mangraph@iD = grapher@rarg, 10D, 8i, 0, 20, 1<D;
Manipulate@mangraph@iD, 8i, 0, 20, 1<DD
**

i

0 p 4

p

2

3 p

4 p

f

2

4

6

8

10

FHr,fL r = 2.00

*legendre.nb ***21**