19.3 The Poisson or Laplace equation in a finite region

19.3.1 Overview

Figure 19.9: Example finite domain $\Omega$ in which the Poisson or Laplace equation is to be solved.

This section will derive the solution of the Poisson equation in a finite region as sketched in figure 19.9. The region will be denoted as $\Omega$, and its boundary by $\delta\Omega$. It will again be assumed that the region is two-di­men­sion­al, leaving the three-di­men­sion­al case to the homework. As shown in figure 19.9, inside the region the Poisson equation applies. In case $f$ $\vphantom0\raisebox{1.5pt}{$=$}$ 0, this becomes the Laplace equation. On the boundary there is some boundary condition that will for now be left arbitrary.

19.3.2 Intro to the solution procedure

The big idea is to relate the finite domain solution to the infinite domain solution derived earlier.

Figure 19.10: Temperature distributions involved in the solution process: (i) $u$ is the desired temperature distribution satisfying the given boundary conditions; (ii) $u_{\rm out}$ is an external temperature distribution whose boundary conditions can be cleverly chosen to achieve various objectives; (iii) $u_{\rm inf}$ is the infinite-domain solution that satisfies no particular boundary conditions on $\delta\Omega$.

For the given heat addition $f$, you can still do the integral $\int{G}f{ \rm d}{V}_{\vec\xi}$ over the domain $\Omega$ to get the infinite domain solution. That solution will now be called $u_{\rm inf}$, as indicated in figure 19.10.

The infinite domain solution satisfies the Poisson equation, but it does not satisfy the boundary conditions. It will turn out that some surface integrals must be added to it to get the boundary conditions right.

The precise form of these integrals can vary. The different versions all give the same, correct, solution $u$ inside the domain $\Omega$. However, they give different answers for the continuation of this solution to outside the domain $\Omega$. So a meaningful discussion of the various possibilities requires consideration of the solution outside the domain. Even though the solution outside domain $\Omega$ is not actually a part of the problem.

The solution outside $\Omega$ will be indicated as $u_{\rm out}$, so the picture becomes as shown in figure 19.10. Since you surely do not want to just make up an arbitrary function $f$ outside $\Omega$, it will be assumed that $f$ $\vphantom0\raisebox{1.5pt}{$=$}$ 0 outside. So $u_{\rm out}$ satisfies the homogeneous Poisson equation, the Laplace equation. And so does the infinite space solution $u_{\rm {inf}}$ outside $\Omega$, for that matter. Only integrating $f$ over the domain $\Omega$ is the same as setting $f$ to zero outside the domain.

19.3.3 Derivation of the integral solution

(book, example 8.2)

The desired integral solution $u$ for the finite-region Poisson solution is a generalization of the infinite domain solution

$$u_{\rm inf}(\vec x) = \int G(\vec x;\vec\xi) f(\vec\xi) { \rm d}V_{\vec\xi},$$

where $G$ is the infinite domain Green’s function derived in the previous section. In the two-di­men­sion­al case discussed here:

$$G = \frac{1}{2\pi}\ln\left\vert\vec x - \vec\xi\right\vert.$$

Since $f$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\nabla^2u$, the expression above represents an important relationship between the infinite-domain solution $u_{\rm inf}$ and the actual, finite-domain, solution $u$:

$$u_{\rm inf}(\vec x) = \int G(\vec x;\vec\xi) \nabla^2_{\vec\xi} u(\vec\xi) { \rm d}V_{\vec\xi}$$

where

$$\nabla^2_{\vec\xi} = \frac{\partial^2}{\partial\xi^2} + \frac{\partial^2}{\partial\eta^2}.$$

is the Laplacian with respect to $\xi$ and $\eta$. George Green discovered that the integral in the right hand side could be simplified into surface integrals using the divergence theorem, and that doing so directly relates $u$ to $u_{\rm {inf}}$.

Figure 19.11: Region of integration of the integral for the infinite-space solution. Note that $\delta\Omega$ is a bounding surface of both dark grey domain $\Omega$ and of the light grey exterior region.

Some caution is needed, however. The Green’s function is infinite when $\vec\xi$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\vec x$, and integrals of infinite functions are not proper. And neither are integrals over infinite regions. You must exclude a very small circle around the point $\vec x$ at which $u_{\rm {inf}}$ is desired from the integration, and also the outside of a very large circle, as indicated in figure 19.11. The correct value for $u_{\rm inf}$ can then be obtained as the limit in which the radius $\varepsilon$ of the small circle becomes zero and the radius $R$ of the large circle becomes infinite. Also, you must think of the integral as really consisting of two separate integrations; one over the dark grey region $\Omega$ and one over the light grey exterior of $\Omega$. The reason is that $u$ and its derivatives are not normally continuous on the surface $\delta\Omega$, and the divergence theorem can only be used for fairly smooth functions.

To simplify the remaining discussion, the origin of the coordinate system will be shifted towards the point $\vec x$ at which $u_{\rm {inf}}$ is desired. The integration coordinate $\vec\xi$ can then be described by polar coordinates $\rho$ and $\varphi$ centered around this point. That simplifies the expression to be evaluated to

$$u_{\rm inf} = \int G \nabla^2_{\vec\xi} u { \rm d}V_{\vec\xi} \qquad\mbox{with}\qquad G = \frac{1}{2\pi}\ln\rho$$

To get a divergence integral, move a $\nabla_{\vec\xi}$ out in front, adding a correction term to make up for the error in doing so:

$$u_{\rm inf} = \int\nabla_{\vec\xi}\left(G\nabla_{\vec\xi... ...i}G\right)\left(\nabla_{\vec\xi}u\right){ \rm d}V_{\vec\xi}.$$

Move a $\nabla_{\vec\xi}$ out in front in the second integral to create another divergence integral, adding another correction term:

$$u_{\rm inf} = \int\nabla_{\vec\xi}\left(G\nabla_{\vec\xi... ... \int\left(\nabla_{\vec\xi}^2G\right)u{ \rm d}V_{\vec\xi}.$$

This final correction term, however, is zero. To see why, remember that the Green’s function $G(\vec x,\vec\xi)$ is the temperature distribution due to a spike of heat at point $\vec\xi$. So $\nabla^2G$ $\vphantom0\raisebox{1.5pt}{$=$}$ 0 everywhere except at the singular point $\vec{x}$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\vec\xi$. And since $\vec\xi$ appears exactly the same way in the Green’s function as $\vec x$, then so is $\nabla_{\vec\xi}^2G$ zero.

The remaining two terms become surface (actually, contour in 2D,) integrals using the divergence theorem. In particular:

$$u_{\rm inf} = - \oint_{\rho=\varepsilon} G \frac{\partial u}{\partial\rho} { \... ... \oint_{\rho=\varepsilon} u \frac{\partial G}{\partial\rho}{ \rm d}S_{\vec\xi}$$

$$+ \oint_{\rho=R} G \frac{\partial u_{\rm out}}{\partial\rho} { \... ...\oint_{\rho=R} u_{\rm out} \frac{\partial G}{\partial\rho} { \rm d}S_{\vec\xi}$$

$$+ \oint_{\delta\Omega} G \left(\frac{\partial u}{\partial n_{\vec... ...{\rm out}\right) \frac{\partial G}{\partial n_{\vec\xi}}{ \rm d}S_{\vec\xi}%$$ (19.9)

To verify this expression, note that the surfaces include the small and big circles, and that $\delta\Omega$ counts as both part of the surface of the dark grey region in figure 19.11 as well as part of the surface of the light grey region. The normal vector $\vec{n}$ on $\delta\Omega$ was taken to stick out of region $\Omega$, which accounts for the additional minus sign in the corresponding $u_{\rm {out}}$ terms. Also $\vec{n}\cdot\nabla_{\vec\xi}$ is according to the total differential of calculus the derivative $\partial/\partial{n}_{\vec\xi}$ in the direction normal to the surface. On the big circle, that is the same as $\partial/\partial\rho$, and on the small circle it is $-\partial/\partial\rho$ since there the outward normal points towards the origin.

The second integral over the small circle is particularly interesting: since $G$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\ln(\rho)/2\pi$, its $\rho$ derivative is $1/2\pi\rho$, which is the inverse of the surface (perimeter) of the circle. So you get

$$\oint_{\rho=\varepsilon} u \frac{\partial G}{\partial\rho}... ..._{\rho=\varepsilon} u \frac{{\rm d}S_{\vec\xi}}{S_{\vec\xi}}.$$

That is just the average of $u$ on the small circle, and it becomes $u$ at the point $\vec x$, (used here as origin,) in the limit that the radius of the small circle becomes zero. So, since this integral simplifies to $u$, all the other integrals in equation (19.9) merely describe the difference between the true solution $u$ and the infinite-domain solution $u_{\rm {inf}}$.

The first integral over the small circle in (19.9) is vanishingly small and can be ignored. To see why, note that it is no larger than the maximum value of the gradient of $u$ on the small circle times

$$\oint_{\rho=\varepsilon} G { \rm d}S_{\vec\xi} = \oint_... ...) \; \rho{ \rm d}\varphi = \ln(\varepsilon) \; \varepsilon$$

and that becomes zero in the limit $\varepsilon\to0$.

There is little that can be done about the integrals over surface $\delta\Omega$. However, the integrals over the big circle in (19.9) still must be evaluated. To do so, you must know something about the behavior of the solution $u_{\rm {out}}$ for large values of $\rho$. In general, it is described by

$$u_{\rm {out}} \sim C_0 \ln\rho + C_1 + \frac{C_2}{\rho} + \ldots$$

In three or more dimensions, the constant $C_0$ is zero. The two integrals over the big circle become, noting that ${\rm d}{S}_{\vec\xi}$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\rho{ \rm d}\varphi$,

$$\oint_{\rho=R} \frac{1}{2\pi} \ln(\rho) \left(\frac{C_0}... ...C_2}{\rho}\right) \frac{1}{2\pi\rho} \;\rho{ \rm d}\varphi$$

which becomes $-C_1$ in the limit $R\to\infty$.

Collecting the results together, the solution for the temperature $u$ at any point $(x,y)$ is:

$$\fbox{$\displaystyle u(\vec x) = \int_{\Omega} G f { ... ...rtial n_{\vec\xi}}\right) { \rm d}S_{\vec\xi} + C_1 $} %$$ (19.10)

where $G$ is the infinite domain Green’s function, $(2\pi)^{-1}\ln\vert\vec x-\vec\xi\vert$, with $\vert\vec x-\vec\xi\vert$ the distance between the point of integration $\vec\xi$ $\vphantom0\raisebox{1.5pt}{$=$}$ $(\xi,\eta)$ and the point $\vec x$ $\vphantom0\raisebox{1.5pt}{$=$}$ $(x,y)$ at which the temperature $u$ is desired. The first integral is therefor the infinite domain solution $u_{\rm inf}$, which has the right values for the added heat $f$, but does not satisfy the correct boundary condition on $\delta\Omega$

19.3.3 Review Questions

1.

Perform the equivalent analysis in the three dimensional case.

Solution pnfd-a

19.3.4 Boundary integral (panel) methods

The previous subsection derived the solution to the Poisson equation in a finite domain. It was given by equation (19.10). This subsection will examine how this solution may be evaluated.

Except for $G(\vec x;\vec\xi)$, all other quantities in the right hand side of equation (19.10) are evaluated at the point of integration $\vec\xi$. For example, $\partial u/\partial n_{\vec\xi}$ stands for the normal derivative $\partial u/\partial n$ evaluated at the boundary point $\xi$ of integration. That means that if you merely know $u$ and the normal derivative $\partial u/\partial n$ on the boundary, you can find $u$ in the interior by taking $u_{\rm out}$ to be zero and doing the integrals above. Unfortunately, a priori at most only one of $u$ (Dirichlet boundary condition) or $\partial u/\partial n$ (Neumann boundary condition) will be known on the boundary.

Various solutions for this problem are possible. A panel method might decide to compute the particular solution where $u_{\rm out}$ is not zero, but has the same values as $u$ on the boundary. The big advantage is then that the second integral in (19.10) drops out, leaving only the last integral as a problem.

A simple panel method will now discretize the boundary in a large number of densely spaced points, and then put a Green’s function at each point. Since each Green’s function corresponds to the addition of a spike of heat at that point, this is called a surface source distribution. The problem remains that the strengths

$$-\left(\partial{u}/\partial{n_{\vec\xi}} -\partial{u}_{\rm {out}}/\partial{n_{\vec\xi}} \right){ \rm d}{S_{\vec\xi}}$$

of these sources are not known, since even if $\partial u/\partial n$ is given on the boundary, $\partial u_{\rm out}/\partial n$ is not. So the strength of each source is an unknown, and an equally large number of equations is needed. These equations can be found from requiring that at that many points, the error in the boundary condition as computed from (19.10) is zero. Put all these equations on a computer and solve. And with the source strength now known, $u$ can then be evaluated at any arbitrary point.

Alternatively, a panel method might decide to compute the solution for the case that $u$ and $u_{\rm out}$ have the same normal derivatives on the boundary. That kills off the source integral, leaving the second integral in (19.10). The quantity $\partial{G}/\partial n_{\vec\xi}$ in this integral is called a dipole. The reason for that name can be understood by writing the definition of the derivative:

$$\frac{\partial G}{\partial n_{\vec\xi}} \equiv \lim_{\De... ...Delta n) - \frac{1}{\Delta n}G(\vec x;\vec\xi) \right). %$$ (19.11)

This shows that a dipole corresponds to an infinitely large source of heat and an infinitely large sink of heat infinitely close together.

19.3.5 Poisson’s integral formulae

The previous subsection showed that the Poisson equation can be solved by using suitable source and/or dipole distributions on the boundary of the domain. However, the strengths of these distributions are not usually known, since they involve both $u$ and its normal derivative on the boundary, and there is only one boundary condition. And if an exterior solution $u_{\rm out}$ is chosen to eliminate one of them, that has the effect of introducing the unknown values of $u_{\rm out}$ or its derivative into the problem. So at least one distribution strength must be found using brute numerical force. Or by brute analytical force, maybe, if the domain is simple.

There is an exception, however, and it occurs for the Dirichlet problem inside a ball (a circle in two dimensions, a sphere in three-dimensions, etcetera.) In that case, suitable distribution strengths can be found by simple means.

The following discussion will restrict itself to the Laplace equation, since the Poisson equation can always be turned into the Laplace equation by subtracting the unbounded space solution $u_{\rm inf}$. This only produces an unimportant change in the inhomogeneous term of the boundary condition. The problem to be solved is then:

$$\nabla^2 u = 0 \mbox{ for } \vert\vec x\vert \le R, \qquad u = g \mbox{ on } \vert\vec x\vert = R,$$

where $g$ is a given function, physically the temperature on the boundary in heat conduction problems, and $R$ is the radius of the ball.

In two dimensions, using polar coordinates, the solution is

$$\fbox{$\displaystyle u(r,\vartheta) = \frac{R^2 - r^2}{2... ...ta} {R^2 - 2 R r\cos(\bar\vartheta-\vartheta) + r^2} $} %$$ (19.12)

and in three dimensions, using spherical coordinates, the solution is

$$\fbox{$\displaystyle \begin{array}{l} \displaystyle ... ...rphi-\varphi)] + r^2 \right\}^{3/2}} \end{array} $} %$$ (19.13)

These results are known as “Poisson’s integral formula” in two, respectively three dimensions.

19.3.6 Derivation

This subsection will derive the two-di­men­sion­al formula above, leaving the three-di­men­sion­al one for the homework. For simplicity, from now on it will be assumed that the ball (i.e. circle in two-dimensions) has unit radius,

$$R = 1$$

It is a simple matter of rescaling $r$ to get back to the formulae for a ball of arbitrary radius.

The integral formula can be derived by a clever selection for the solution $u_{\rm {out}}$ outside the circle in the integral solution (19.10). In particular, the trick is to take

$$u_{\rm out}(r,\vartheta) = A u(\bar r,\vartheta) \qquad\mbox{where}\qquad \bar r = \frac{1}{r} %$$ (19.14)

Here $A$ is a constant still to be selected. Note that if $\bar{r}$ $\raisebox{-.3pt}{$\leqslant$}$ 1 then $r$ $\raisebox{-.5pt}{$\geqslant$}$ 1: these rules turn solutions inside the ball into solutions outside the ball. The transformation $\bar{r}$ $\vphantom0\raisebox{1.5pt}{$=$}$ $1/r$ is called an inversion with respect to the surface of the unit ball.

The first thing to show is that $u_{\rm {out}}$ satisfies the Laplace equation. The integral solution (19.10) does not apply otherwise. The Laplacian,

$$\nabla^2 u_{\rm out} = \frac{\partial^2 u_{\rm out}}{\pa... ...frac{1}{r^2} \frac{\partial^2 u_{\rm out}}{\partial \theta^2}$$

must be zero.

To show that this is so, first differentiate (19.14) once, using the chain rule to convert the $\bar{r}$ derivatives of $u$ to $r$ derivatives:

$$\frac{\partial u_{\rm out}}{\partial r} = A \frac{\parti... ...al r}= - A \frac{\partial u}{\partial \bar r} \frac{1}{r^2}$$

Differentiate this once more to get the second derivative. Note that you now have to use the product rule of differentiation to differentiate the factors. And you need again the chain rule for differentiating the first factor. You get

$$\frac{\partial^2 u_{\rm out}}{\partial r^2} = A \frac{\p... ...}{r^2} + A \frac{\partial u}{\partial \bar r} \frac{2}{r^3}$$

Also,

$$\frac{\partial^2 u_{\rm out}}{\partial \theta^2} = A \frac{\partial^2 u}{\partial \theta^2}$$

If you plug these derivatives into the Laplacian given above, you get

$$\nabla^2 u_{\rm out} = A \frac{1}{r^4} \left[ \frac{\p... ...r r} + r^2 \frac{\partial^2 u}{\partial \theta^2} \right]$$

Since $r$ $\vphantom0\raisebox{1.5pt}{$=$}$ $1/\bar{r}$, you recognize the Laplacian of $u$ inside the square brackets. That is zero because $u$ satisfies the Laplace equation. Then you see that so does $u_{\rm {out}}$.

Now the idea is to try to choose the constant $A$ so that the integral solution (19.10) only involves the given values of $u$ on the boundary. In particular, the normal derivative of $u-u_{\rm {out}}$ must be eliminated. Now for a spherical boundary, the normal derivative is the radial derivative. And on the surface of the ball, $r$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\bar{r}$ $\vphantom0\raisebox{1.5pt}{$=$}$ 1. So on the boundary, using the expressions above,

$$u_{\rm out}(1,\theta) = A u(1,\theta) \qquad \frac{\pa... ...eta) = -A \frac{\partial u_{\rm out}}{\partial r}(1,\theta)$$

Note that in the final term, the first independent variable in $u$ has been renamed simply $r$. It does not make a difference what you call the independent variable of a function; we just used a bar on it when we were treating $u$ at one location to define $u_{\rm {out}}$ at another location. The bar was merely to keep the two locations apart.

For $\partial(u-u_{\rm {out}})/\partial{r}$ to vanish on the surface of the sphere. according to the above equations you need to take $A$ $\vphantom0\raisebox{1.5pt}{$=$}$ $-1$. In that case, $u-u_{\rm {out}}$ on the sphere equals $2u$, and $u$ is the given function $g$ on the surface of the sphere. So the integral solution (19.10) becomes

$$u(\vec x) = \oint_{\vert\vec\xi\vert=1} 2 g \frac{\parti... ... G = \frac{1}{2\pi}\ln\left\vert\vec x-\vec\xi\right\vert. %$$ (19.15)

The above solution $u$ is completely in terms of the given function $g$. So the Dirichlet problem has been solved.

But of course you want to clean it up. You would like the solution of a problem in a circle to be in terms of polar coordinates. So set

$$\vec x = r {\hat\imath}_r \mbox{ with } {\hat\imath}_r= ... ...t(\begin{array}{c}\cos\varphi \sin\varphi\end{array}\right)$$

for the point $\vec x$ at which the temperature $u$ is desired and the point of integration $\vec\xi$ respectively. Then the surface element ${\rm d}{S_{\vec\xi}}$ in the integral over the circle perimeter is $\rho{ \rm d}\varphi$, and $\rho$ $\vphantom0\raisebox{1.5pt}{$=$}$ 1 on the circle.

Also, the derivative $\partial{G}/\partial{n}_{\vec\xi}$ normal to the circle is simply $\partial{G}/\partial\rho$. G is a function of the distance $d$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\vert\vec{x}-\vec\xi\vert$ between the points $\vec{x}$ and $\vec\xi$; in particular $G$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\ln(d)/2\pi$ in two dimensions. You can write

$$d^2 = (r{\hat\imath}_r-\rho{\hat\imath}_\rho)\cdot(r{\hat\... ... r^2 - 2r\rho{\hat\imath}_r\cdot{\hat\imath}_\rho + \rho^2, %$$ (19.16)

whose derivative with respect to $\rho$ equals

$$\frac{{\rm d}d^2}{{\rm d}\rho} = \frac{- 2r\rho{\hat\imath}_r\cdot{\hat\imath}_\rho + 2\rho^2}{\rho}$$

or getting rid of the ugly dot product term using the expression (19.16) for $d^2$,

$$\frac{{\rm d}d^2}{{\rm d}\rho} = \frac{\rho^2 - r^2 + d^2}{\rho}$$

So you can write using the chain rule that

$$\frac{\partial G}{\partial n_{\vec\xi}} = \frac{{\rm d... ... \frac{{\rm d}G}{{\rm d}d} \frac{\rho^2 - r^2 + d^2}{2d\rho}$$

Plug in the expression $G$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\ln(d)/2\pi$ for the two-di­men­sion­al Green’s function, and note that $\rho$ $\vphantom0\raisebox{1.5pt}{$=$}$ 1 on the circle to get:

$$\frac{\partial G}{\partial n_{\vec\xi}} = \frac{1 - r^2}{4\pi d^2} + \frac{1}{4\pi}$$

Plug that into the integral expression (19.15) for $u(\vec{x})$, taking $d^2$ from (19.16) with ${\hat\imath}_r\cdot{\hat\imath}_\rho$ equal to $\cos(\varphi-\vartheta)$, to get

$$u(r,\vartheta) = \frac{1 - r^2}{2\pi} \oint \frac{g(\var... ...2} + \frac{1}{2\pi} \oint g(\varphi){ \rm d}\varphi + C_1.$$

The two final terms are just constants, and they cancel each other. The reason is that

$$C_1=u_{\rm out}(\infty,\vartheta)=-u(0,\vartheta)$$

The mean value theorem, proved in {D.3}, says that $u(0,\vartheta)$ equals the average of $g$ $\vphantom0\raisebox{1.5pt}{$=$}$ $u$ on the circle.

Also, to allow for the case that the radial coordinate $r$ is not normalized with the circle radius $R$, you want to replace $r$ in the above result with $r/R$. That produces the Poisson integral as stated in the previous subsection.

19.3.6 Review Questions

1.

Find a suitable solution $u_{\rm {out}}$ outside the sphere in three dimensions. Show that it satisfies the Laplace equation.

Solution pnifd-a

2.

Derive the Poisson integral formula in three dimensions as given in the previous subsection.

Solution pnifd-b

19.3.7 The integral formula for the Neumann problem

The Neumann problem in two dimensions is:

$$\nabla^2 u = 0 \mbox{ for } r \le R, \qquad \frac{\partial u}{\partial r} = g \mbox{ on } r = R,$$

This corresponds physically to a problem where the heat flux instead of the temperature is described on the boundary. The solution is

$$\fbox{$\displaystyle u(r,\vartheta) = \frac{-1}{2\pi} ... ...\varphi-\vartheta) + r^2}{R^2}{ \rm d}\varphi + C_1 $} %$$ (19.17)

Note that there is only a proper solution for $u$ if

$$\oint g(\varphi){ \rm d}\varphi = 0$$

If you put in an invalid $g$, you will get a $u$, but it will not have heat flux $g$ through the boundary. In particular, putting in a nonzero constant for $g$ will produce $u$ $\vphantom0\raisebox{1.5pt}{$=$}$ 0 and no heat flux. It can be seen from the above expression that the undetermined constant $C_1$ is the temperature at the center of the circle.

The derivation of the formula above is similar to the one in the previous subsection. You will be disappointed to learn that you must miss doing it in the homework. The same story does not work in three dimensions since you cannot get rid of the unknown surface values of $u$ in both the source and dipole distributions. In two dimensions, however, if you take $u_{\rm out}(r,\vartheta)$ $\vphantom0\raisebox{1.5pt}{$=$}$ $u(1/r,\vartheta)$, the dipole strength is zero and only the source integral remains:

$$u(\vec x) = - \oint_{\vert\vec\xi\vert=1} 2 g G{ \rm d}S_{\vec\xi} + C_1 \mbox{ with } G = \frac{1}{2\pi}\ln d$$

and use of the expression (19.16) for $d$ gives the stated result.

19.3.8 Smoothness of the solution

One important qualitative conclusion that can be drawn from the various results of the previous subsections is that the solution of a Laplace equation problem is infinitely smooth in the interior of the region $\Omega$ in which it applies.

For example, consider the derived expression for $u$ if the exterior solution is zero:

$$u(\vec x) = \oint_{\delta\Omega} u \frac{\partial G}{\pa... ...G \frac{\partial u}{\partial n_{\vec\xi}}{ \rm d}S_{\vec\xi}$$ (19.18)

If you take derivatives of $u$ with respect to the components of $\vec{x}$, you will be differentiating $G(\vec x;\vec\xi)$ inside the integral. And $G$ has infinitely many finite derivatives away from the singular point $\vec\xi$ $\vphantom0\raisebox{1.5pt}{$=$}$ $\vec{x}$, in other words, away from the boundary.

So, if $u$ and $\partial{u}/\partial{n}$ are merely integrable on the boundary, which still allows them to be quite singular, the solution at every point in the interior will have infinitely many continuous derivatives.

It is somewhat different for the Poisson equation, since if the forcing $f$ has a singularity at some point, then so will the solution $u$. But still the solution for $u$ will be less singular than $f$ is. For example, in two-dimensions a delta function in $f$, whose square is not integrable, produces a logarithmic Green’s function, for which every power is integrable over the singular point. In general, it can be seen from Fourier solution of the Poisson problem that $u$ will in general have two more square integrable derivatives than $f$. (Assuming that lack of decay of $u$ at large distances is not a factor or subtracted out first.)