Preface
This section concerns about Laplace's equation in spherical coordinates.
Return to computing page for the first course APMA0330
Return to computing page for the second course APMA0340
Return to Mathematica tutorial for the first course APMA0330
Return to Mathematica tutorial for the second course APMA0340
Return to the main page for the first course APMA0330
Return to the main page for the second course APMA0340
Return to Part VI of the course APMA0340
Introduction to Linear Algebra with Mathematica
Glossary
The Laplace equation for a function \( u (r, \phi , \theta ) \) is given by
In the n-dimensional space, for n ≥ 3, the uniqueness of the exterior problem does not valid even in the class of smooth functions bounded outside the given closed domain.
Let Ω be a unit sphere in the n-dimensional Eucledian space ℝn, defined by the equation
\begin{equation} \label{EqSphere.1}
\nabla^2 u (r, \phi , \theta ) = \frac{1}{r^2} \,\frac{\partial}{\partial r} \left( r^2 \frac{\partial u}{\partial r} \right) + \frac{1}{r^2 \sin\phi} \, \frac{\partial}{\partial \phi} \left( \sin \phi \, \frac{\partial u}{\partial \phi} \right) + \frac{1}{r^2 \sin^2 \phi} \, \frac{\partial^2 u}{\partial \theta^2} =0 .
\end{equation}
When considering only axisymmetric solutions, that is, solutions which depend on r and φ but not on θ, the Laplace equation is reduced to
\begin{equation} \label{EqSphere.2}
\nabla^2 u (r, \phi ) = \frac{1}{r^2} \,\frac{\partial}{\partial r} \left( r^2 \frac{\partial u}{\partial r} \right) + \frac{1}{r^2 \sin\phi} \, \frac{\partial}{\partial \phi} \left( \sin \phi \, \frac{\partial u}{\partial \phi} \right) =0 .
\end{equation}
In particular, the electrostatic potential satisfies the above equation.
Example: Consider the outer Dirichlet's problem
\[
\Delta u({\bf x}) = 0, \quad {\bf x} \in \Omega = \left\{ {\bf x}\in \mathbb{R}^n \,:\, r = \| {\bf x}\| > 1 \right\} , \qquad u = 1 \quad\mbox{on the boundary}\quad \partial\Omega .
\]
Here \( \displaystyle \| {\bf x} \|^2 = x_1^2 + x_2^2 + \cdots + x_n^2 \) is the square of the the Euclidean norm.
In the space of continuous, but unbounded functions, the given boundary value problem has infinite many solutions.
Indeed, the following one-parameter family of functions satisfies all conditions of the given Dirichlet's problem:
\[
v_{\lambda} (r) = \begin{cases}
\lambda + \frac{1-\lambda}{r^{n-2}} , & \ n > 2, \\
1 + \left( 1 - \lambda \right) \ln r , & \ n = 2 ,
\end{cases}
\]
where −∞ < λ < ∞. Note that eqch function vλ is neither bounded nor tends to zero at infinity.
■ Example: Consider teh harmonic function
\[
u = u(r) = r^{2-n} , \qquad n > 2.
\]
This bounded function is a solution of the homogebeous Dirichlet's problem:
\[
\Delta u = \frac{\partial^2 u}{\partial x_1^2} + \frac{\partial^2 u}{\partial x_2^2} + \cdots + \frac{\partial^2 u}{\partial x_n^2} =0 \quad \mbox{outside} \quad r> 1, \qquad u = 1 \quad \mbox{on the boundary} \quad \| {\bf x} \| = 1.
\]
■
\begin{equation} \label{EqSphere.3}
\| {\bf x} \| = 1 , \qquad\mbox{having area} \quad \omega_n = \frac{2\,\pi^{n/2}}{\Gamma \left( \frac{n}{2} \right)} ,
\end{equation}
where Γ(z) is the gamma function of Euler.
If u = u(x), x = (x1, x2, … , xn) is a regular harmonic function in a domain E of x-space with boundary ∂E, then the function
The harmonic function \eqref{EqSphere.4} is called regular in the exterior region E*. That is, we define regularity in a domain E including the infinity as follows: we invert E with respect to sphere
\( \| {\bf x} \| = 1 , \) transfering it into a bounded domain G*.
\begin{equation} \label{EqSphere.4}
v ({\bf x}) = \frac{1}{r^{n-2}} \left( \frac{x_1}{r^2} , \frac{x_2}{r^2} , \ldots , \frac{x_n}{r^2} , \right)
\end{equation}
also satisfies the Laplace equation Δu = 0 and is regular in the region E* obtained from E by inversion with respect to the unit sphere \( \| {\bf x} \| = 1 . \)
The map x ↦ x*, where
Note that if x ≠ {0, {∞}, then x* lies on the ray from the origin determined by x.
\begin{equation} \label{EqSphere.5}
{\bf x}^{\ast} = \begin{cases} \frac{\bf x}{\| {\bf x} \|^2} , & \ \mbox{if } \ {\bf x} \ne 0, \infty , \\
0 , & \ \mbox{if } \ {\bf x} = \infty , \\
\infty , & \ \mbox{if } \ {\bf x} = 0, \end{cases}
\end{equation}
is called the inversion of ℝn ∪ {∞} relative to the unit sphere.
Temperatures in a Sphere
The steady temperature distribution u(x) inside the sphere r = a, in spherical polar coordinates, satisfies \( \nabla^2 u =0 . \) If we heat the surface of the sphere so that \( u = f(\theta ) \) on r = a for some given function \( f(\theta ) , \) what is the temperature distribution within the sphere?
The equation and boundary conditions do not depend on φ, so we know that u is of the form
\[
u(r, \theta ) = \sum_{n\ge 0} \left[ A_n r^n + \frac{B_n}{r^{n+1}} \right] P_n (\cos \theta )
\]
for some constants An and Bn. Furthermore we expect u to be finite at r = 0 so that \( B_n =0 . \) We find the coefficients An by evaluating the sum at r =a:
\[
f( \theta ) = \sum_{n\ge 0} A_n a^n \, P_n (\cos \theta ) .
\]
We can find An using the orthogonality of the Legendre polynomials. However, the integration is done with respect to x but not \( \cos ( \theta ) . \) Setting \( x= \cos ( \theta ) , \) we get \( {\text d}x= - \sin ( \theta ) \,{\text d} \theta . \) The interval of integration [-1,1] becomes \( [-\pi , 0 ] . \) Multiply through by
\( -\sin ( \theta ) \,P_m (\cos \theta ) \) and integrate in θ to obtain
\begin{align*}
\int_{-\pi}^0 -\sin ( \theta ) \,P_m (\cos \theta ) \,f(\theta )\, {\text d} \theta &= \sum_{n\ge 0} \left( a^n A_n \right) \int_{-\pi}^0 -\sin ( \theta ) \,P_m (\cos \theta ) \,P_n (\theta )\, {\text d} \theta
\\
&= \sum_{n\ge 0} \left( a^n A_n \right) \int_{-1}^1 P_n (x)\, P_m (x) \, {\text d}x = \frac{2\,a^m A_m}{2m+1} .
\end{align*}
So
\[
u(r, \theta ) = \sum_{n\ge 0} \left( n + \frac{1}{2} \right) \left( \frac{r}{a} \right)^n \, P_n (\cos \theta ) \, \int_0^{\pi} f(\nu )\,P_n (\cos \nu )\,\sin \nu \,{\text d}\nu .
\]
Example: Consider Laplace's equation exterior to a sphere of radius a, subject to some boundary condition on the sphere. The full problem statement is given below
\[
\frac{1}{r^2} \,\frac{\partial}{\partial r} \left( r^2 \frac{\partial \Phi}{\partial r} \right) + \frac{1}{r^2 \sin\phi} \, \frac{\partial}{\partial \phi} \left( \sin \phi \, \frac{\partial \Phi}{\partial \phi} \right) =0 , \quad r>a \quad\mbox{and} \quad 0 \le \phi \le \pi ,
\]
\[
\Phi (a,\phi ) = \left\{
\begin{array}{ll}
V_0 , & \ \mbox{for $\phi \in \left( 0 , \frac{\pi}{2} \right)$} \\
-V_0 , & \ \mbox{for $\phi \in \left( \frac{\pi}{2} , \pi \right)$. }
\end{array}
\right.
\]
\[
\Phi (r, \phi ) \, \mapsto \, 0 \quad \mbox{when} \quad r \mapsto 0.
\]
Because of the condition at \( \infty , \) only the negative powers of r may be used in this region exterior to a sphere. We superpose all of those solutions to get
\[
\Phi (r, \phi ) = \sum_{n \ge 0} B_n r^{-n-1} \,P_n \left( \cos \phi \right) .
\]
Then imposing the boundary condition, we get
\[
\Phi (a, \phi ) = V_0 \,\mbox{sign}\left( \phi - \frac{\pi}{2} \right) = \sum_{n \ge 0} B_n a^{-n-1} \,P_n \left( \cos \phi \right) .
\]
Previously, we found coefficients of Legendre expansion for the signum function, so we get
\[
B_n a^{-n-1} = V_0 \left[ P_{n-1} (0) - P_{n+1} (0) \right] , \qquad n=1,2 \ldots .
\]
Example: As our final example, we consider the region between two concentric spheres, with radii \( a \mbox{ and } b , \quad 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.
\[
\frac{1}{r^2} \,\frac{\partial}{\partial r} \left( r^2 \frac{\partial \Phi}{\partial r} \right) + \frac{1}{r^2 \sin\phi} \, \frac{\partial}{\partial \phi} \left( \sin \phi \, \frac{\partial \Phi}{\partial \phi} \right) =0 , \quad r \in (a, b) \quad\mbox{and} \quad 0 \le \phi \le \pi ,
\]
\[
\Phi (a, \phi ) = g(\phi ), \qquad \Phi (b, \phi ) = h(\phi ) .
\]
We begin by expanding both g and h in Legendre polynomials. That will make our task easier
\[
g (\phi ) = \sum_{n\ge 0} c_n P_n \left( \cos \,\phi \right) , \qquad c_n = \frac{2n+1}{2} \, \int_0^{\pi} g(\phi )\, P_n \left( \cos \,\phi \right) \sin \,\phi \,{\text d} \phi ;
\]
\[
h (\phi ) = \sum_{n\ge 0} d_n P_n \left( \cos \,\phi \right) , \qquad d_n = \frac{2n+1}{2} \, \int_0^{\pi} h(\phi )\, P_n \left( \cos \,\phi \right) \sin \,\phi \,{\text d} \phi .
\]
The coefficients cn and dn are known from the known boundary functions g and h. The separated solutions are \( r^n P_n \left( \cos\,\phi \right) \) and \( r^{-n-1} P_n \left( \cos\,\phi \right) . \) The domain of the present problem, \( r \in (a, 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 Φ is then obtained by superposition:
\[
\Phi (r, \phi ) = \sum_{n\ge 0} \left( a_n r^n + b_n r^{-n-1} \right) P_n \left( \cos \,\phi \right) .
\]
We now impose the two boundary conditions, so
\[
\Phi (a, \phi ) = \sum_{n\ge 0} \left( a_n a^n + b_n a^{-n-1} \right) P_n \left( \cos \,\phi \right) = g(\phi ) = \sum_{n\ge 0} c_n P_n \left( \cos \,\phi \right) ,
\]
and
\[
\Phi (b, \phi ) = \sum_{n\ge 0} \left( a_n b^n + b_n b^{-n-1} \right) P_n \left( \cos \,\phi \right) = h(\phi ) = \sum_{n\ge 0} d_n P_n \left( \cos \,\phi \right) ,
\]
By equating the coefficients of corresponding terms in the Legendre expansions, we get
\[
a_n a^n + b_n a^{-n-1} = c_n , \qquad \mbox{and} \qquad
a_n b^n + b_n b^{-n-1} = d_n , \qquad n=0,1,2,\ldots .
\]
These constitute two linear algebraic equations for each pair an and bn, which we solve with the aid of Mathematica:
\[
a_n = \frac{b^{n+1} d_n - a^{n+1} c_n}{b^{2n+1} - a^{2n+1}} \qquad\mbox{and}\qquad
b_n = \frac{\left( ab\right)^{n+1} \left( b^n c_n - a^{n} d_n \right)}{b^{2n+1} - a^{2n+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:
\[
g(\phi )= H(\phi - \pi /2) , \qquad h(\phi ) \equiv 0,
\]
where H is the Heaviside function. Now we calculate coefficients:
coeff[n_] := ((n ( 2 n + 3) LegendreP[n-2, 0] - (2 n +1) LegendreP[n, 0]) - (n+1) (2 n -1) LegendreP[n+2, 0])/(2 (2 n -1) (2 n +3))
We put first 30 coefficients, starting at n = 1 into array coef2:
Module[{i}, coef2={0.5};
Do[coef2 = Append[coef2, N[coeff[i]]], {i, 2, 30}]]
The module below defines the kth partial sum of the Legendre expansion and assigns it to legsum2:
legsum2[x_, k_]:= Module[{i}, 0.25+Sum[N[coef2[[i]]]*LegendreP[i, x],{i, 1, k}]]
So \( c_n = coef2[[n]] , \quad d_n =0 . \)
a[n_] := -a^(n+1) coeff[n]/(b^(2*n+1) - a^(2*n+1))
b[n_] := a^(n+1)*b^(2*n+1) coeff[n]/(b^(2*n+1) - a^(2*n+1))
We now define the kth partial sum of the solution Φ:
b[n_] := a^(n+1)*b^(2*n+1) coeff[n]/(b^(2*n+1) - a^(2*n+1))
Phisum3[r_,phi_,k_] := Sum[(a[i]*r^(i) + b[i]*r^(-i-1) )*LegendreP[i, Cos[phi]], {i,0,k} ]
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; b = 5;
Next we define a function grapher[r,k] which uses the kth partial sum to construct a plot of potential versus φ on the sphere of radius r:
grapher[r_,k_] := Plot[ Phisum3[r,phi,k] , {phi,0,Pi}], AxesLabel-> {"\[Phi]", "\[CapitalPhi][r,\[Phi]]"}, PlotLabel-> Row[{"r= " , PaddedForm[N[r], {5,2}]}], PlotRange->{0, 1.1}, Ticks -> {{0, Pi/4, Pi/2 , 3*Pi/4 , Pi}, Automatic}]
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.
GraphicsGrid[Table[{grapher[a+i*(b-a)/5, 10], grapher[a+(i+1)*(b-a)/5, 10]} , {i,0,5,2}]]
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.
DynamicModule[{mangraph, i, rarg},
Do[rarg = a + ((b-a)/20)*i; mangraph[i] = grapher[rarg, 10], {i,0,20,1}];
Manupulate[mangraph[i], {i,0,20,1}]]
Do[rarg = a + ((b-a)/20)*i; mangraph[i] = grapher[rarg, 10], {i,0,20,1}];
Manupulate[mangraph[i], {i,0,20,1}]]
- Axler, S., Bourdon, P., Harmonic Function Theory, Second edition, 2000.
Return to Mathematica page
Return to the main page (APMA0340)
Return to the Part 1 Matrix Algebra
Return to the Part 2 Linear Systems of Ordinary Differential Equations
Return to the Part 3 Non-linear Systems of Ordinary Differential Equations
Return to the Part 4 Numerical Methods
Return to the Part 5 Fourier Series
Return to the Part 6 Partial Differential Equations
Return to the Part 7 Special Functions