# Laplace PDE inside a disk

I can't get Mathematica to solve this standard textbook PDE, which is Laplace inside a disk of some radius. One of the boundary conditions needed is that the solution is finite (bounded) in center of disk, and I do not know how specify this boundary condition. I do not know if this is why Mathematica can't solve it. Most examples in help are for NDSolve and finite elements, and I am looking to analytical solution.

This is the formal description of the PDE

Solve for $u\left( r,\theta\right)$ the PDE

\begin{align*} \frac{\partial^{2}u}{\partial r^{2}}+\frac{1}{r}\frac{\partial u}{\partial r}+\frac{1}{r^{2}}\frac{\partial^{2}u}{\partial\theta^{2}} & =0\\ 0 & \leq r\leq a\\ 0 & <\theta\leq2\pi \end{align*}

Boundary conditions

\begin{align*} u\left( a,\theta\right) & =f\left( \theta\right) \\ \left\vert u\left( 0,\theta\right) \right\vert & <\infty\\ u\left( r,0\right) & =u\left( r,2\pi\right) \\ \frac{\partial u}{\partial\theta}\left( r,0\right) & =\frac{\partial u}{\partial\theta}\left( r,2\pi\right) \end{align*}

The standard solution using seperation of variables by hand gives

\begin{align*} u\left( r,\theta\right) & =A_{0}+\sum_{n=1}^{\infty}r^{n}\left( c_{n}\cos n\theta+b_{n}\sin n\theta\right) \\ A_{0} & =\frac{1}{2\pi}\int_{0}^{2\pi}f\left( \theta\right) d\theta\\ c_{n} & =\frac{1}{\pi a^{n}}\int_{0}^{2\pi}f\left( \theta\right) \cos\left( n\theta\right) d\theta\qquad n>0\\ b_{n} & =\frac{1}{\pi a^{n}}\int_{0}^{2\pi}f\left( \theta\right) \sin\left( n\theta\right) d\theta\qquad n>0 \end{align*}

This is what I tried

ClearAll[u,theta,r,a];
pde=D[u[r,theta],{r,2}]+1/r D[u[r,theta],r]+1/r^2 D[u[r,theta],{theta,2}]==0;
bc=u[a,theta]==f[theta];
sol=DSolve[{pde,bc},u[r,theta],{r,theta},Assumptions->a<r&&a>0&& 0<theta<=2 Pi]


Now I tried adding the bounded condition, hoping it will help. But I do not know what the syntax should be

ClearAll[u,theta,r,a];
pde=D[u[r,theta],{r,2}]+1/r D[u[r,theta],r]+1/r^2 D[u[r,theta],{theta,2}]==0;
bc={u[a,theta]==f[theta],Abs[u[0,theta]]<Infinity}
sol=DSolve[{pde,bc},u[r,theta],{r,theta},Assumptions->a<r&&a>0&& 0<theta<=2 Pi]


Ok, so I tried now to add the periodic conditions. But this did not help

ClearAll[u,theta,r,a];
pde=D[u[r,theta],{r,2}]+1/r D[u[r,theta],r]+1/r^2 D[u[r,theta],{theta,2}]==0;

bc={u[a,theta]==f[theta],
u[r,0]==u[r,2 Pi],
Derivative[0,1][u][r,0]==Derivative[0,1][u][r,2 Pi]}

sol=DSolve[{pde,bc},u[r,theta],{r,theta},Assumptions->a<r&&a>0&& 0<theta<=2 Pi]


Question is: Is Mathematica not able to solve it, because it is missing the bounded conditions part? If so, how to specify this? Any one could get Mathematica to solve this analytically?

This is a very standard PDE problem and its analytical solution is all over the net. The above solution I show is the standard text book solution.

11.3 on windows 7

• The Dirichlet problem under consideration can be solved in Mathematica by making use of Poisson kernel (see en.wikipedia.org/wiki/Poisson_kernel ). – user64494 Mar 16 '18 at 9:55
• @user, maybe write an answer demonstrating how to use the Poisson kernel, should you find time? – J. M.'s technical difficulties Mar 17 '18 at 11:26

Once again, since DSolve doesn't work well at the moment, I'd like to post an answer based on finite Fourier transform:

(* Definition of finiteFourierTransform etc. are not included in this post,
pde = D[u[r, theta], {r, 2}] + 1/r D[u[r, theta], r] +
1/r^2 D[u[r, theta], {theta, 2}] == 0;

bc = {u[a, theta] == f[theta], u[r, 0] == u[r, 2 Pi],
Derivative[0, 1][u][r, 0] == Derivative[0, 1][u][r, 2 Pi]};

Format@finiteFourierTransform[f_, __] := ℱ[f]
Format@theta := θ

finiteFourierTransform[{pde, bc[[1]]}, {theta, 0, 2 Pi}, n]


We substitute the periodic b.c.s in, and solve the resulting system:

% /. Rule @@@ bc[[2 ;; 3]]

tset = % /. HoldPattern@finiteFourierTransform[f_ /; ! FreeQ[f, u], __] :> f
tsol = DSolve[tset, u[r, theta], r][[1, 1, -1]]


One constant C[1] remains, but we still have the bounded condition / finiteness condition at hand. With a bit of transforming:

tsolcollect = Collect[Simplify[tsol // TrigToExp, n ∈ Integers], r^_]


It becomes clear that the coefficient of r^-n should be 0, because $$\lim_{r\to 0^+} r^{-n}=\infty$$ for $$n>0$$:

tsolfinal = tsolcollect /. Flatten@Solve[tsolcollect[[1]] == 0, C[1]] // Simplify


The last step is to transform back:

sol = inverseFiniteFourierTransform[tsolfinal, n, {theta, 0, 2 Pi}, Re]
sol // transformToIntegrate


Though it's troublesome to make Mathematica simplify sol further, it's not hard to notice this answer is equivalent to the one in your question.

sol can be used for further calculation of course. For example, when $$f=\sin(3\theta)$$:

ReleaseHold[sol /. f -> (Sin[3 #] &) /. C -> 7 /. a -> 2] // ComplexExpand // Simplify
(* 1/8 r^3 Sin[3 theta] *)
RevolutionPlot3D[%, {r, 0, 2}, {theta, 0, 2 Pi}, PlotRange -> All]


• Sorry, but C -> 7 (as well as several other steps) is built on the sand. – user64494 Mar 16 '18 at 21:02
• @user64494 Er… what do you mean by "build on the sand"? Taking the first several terms of a series to approximate its value is a common practice. (f -> (Sin[3 #] & happens to be a special case that only finiteFourierTransform[Sin[3 theta], {theta, 0, 2 Pi}, 3] is non-zero so the final solution is not a series though. ) – xzczd Mar 17 '18 at 3:44

This is very strange. I suspect Mathematica desktop is being updated/patched from the cloud when one is connected to the internet, because I just tried again the first attempt I showed in my question above and now DSolve return the solution to the Laplace PDE !

I do have internet setting in preference to allow Wolfram access to do updates. But I thought this is for documentations and things like that only. I do not understand how this could have happened. But Mathematica 11.3 can actually solve this PDE.

code again (same as in my question)

ClearAll[u,theta,r,a];
pde=D[u[r,theta],{r,2}]+1/r D[u[r,theta],r]+1/r^2 D[u[r,theta],{theta,2}]==0;
bc=u[a,theta]==f[theta];
sol=DSolve[{pde,bc},u[r,theta],{r,theta},Assumptions->a<r&&a>0&&0<theta<=2 Pi]
sol=sol/.K[1]->n


The third attempt I had in my question, still does not work, but the first works.

I assume someone tried that also and it did not work for them, else they would have screamed at me saying that it works. Very strange. But happy that 11.3 can solve this.