# How to solve a PDE with periodic boundary conditions in one of the variables?

I need to solve a PDE where one of the variables is an angle, so I need to know how to deal with periodic boundary conditions.

As a warm up, I am trying to solve the Helmholtz equation in polar coordinates, where the solutions should be Bessel functions multiplied by complex exponentials.

I tried the following code:

 NDSolve[{x D[u[x, t], x] + x^2 D[u[x, t], x, x] + D[u[x, t], t, t] + x^2 u[x, t] == 0, u[b, t] == 1, (D [u[x, t], x] /. x -> b) == 0.5, u[x, 0] == u[x, 2 Pi], (D[u[x, t], t] /. t -> 0) == (D[u[x, t], t] /. t -> 2 Pi)}, u[x, t], {x, b, 100}, {t, 0, 2*Pi}]


where i have set b=0.0001 because setting b=0 gives rise to singularities for some reason. One solution to this equation is BesselJ[0,x] + BesselJ[1,x], but it is not the only one: adding any combination of higher order Bessel functions one gets infinite other solutions.

The weird thing is that I noticed, by plotting the result, that Mathematica is giving me only BesselJ[0,x] as the solution! Anyone knows what might be happening?

• u[b, t] == 1 guarantees that there is no azimuthal dependence, which in turn allows only BesselJ[0, x] as a solution. Also, Bessel's equation is singular at x == 0, and b > 0 is necessary to avoid it. – bbgodfrey Jun 6 '18 at 1:30
• I'm voting to close this question as off-topic because the issue it raises is not really a Mathematica issue but a matter of the OP not having grasped the mathematics involved. – m_goldberg Jun 6 '18 at 19:59
• @bbgodfrey Would you consider turning your comment into a short answer? I agree with the voters, that this is basically not a question about Mathematica, but a misunderstanding of the PDE and numerical solving in general. On the other hand, the Helmoltz equation is not that uncommon and we answered many strongly math-related questions in the past. Therefore, I'm a bit hesitant to cast the final close vote. – halirutan Jun 8 '18 at 0:30
• @halirutan I would be happy to do so in a few days, when I have some time. Congratulations on your election as a Moderator. – bbgodfrey Jun 8 '18 at 1:25
• @bbgodfrey Thanks and thank you for considering an answer :) – halirutan Jun 8 '18 at 1:36

As I remarked in comments above, the code in the question actually represents the Laplace equation, not the Helmholtz equation. It is necessary to set b > 0, because the equation is singular at x == 0, as can be seen by solving for D[u[x, t], {x, 2}]. Why it should yield only the zero-order Bessel function can be seen as follows. Let u[x, t] be represented by f[x] Cos[n t]], with n an arbitrary integer. Inserting this into the PDE gives

Expand[#/(Cos[n t] x^2) & /@
FullSimplify[
(x D[u[x, t], x] + x^2 D[u[x, t], x, x] + D[u[x, t], t, t] + x^2 u[x, t] == 0)
/. u -> Function[{x, t}, f[x] Cos[n t]]]
]

(* f[x] - (n^2 f[x])/x^2 + f'[x]/x + f''[x] == 0 *)


which is Bessel's equation of order n. Because the x == b boundary conditions are independent of t, only the n == 0 solution is permitted. That the PDE cannot have as an answer BesselJ[0,x] + BesselJ[1,x] also can be demonstrated by.

FullSimplify[
(x D[u[x, t], x] + x^2 D[u[x, t], x, x] + D[u[x, t], t, t] + x^2 u[x, t] == 0)
/. u -> Function[{x, t}, BesselJ[0, x] + BesselJ[1, x]]]

(* BesselJ[1, x] == 0 *)


Incidentally, numerically solving Laplace's equation as an initial value problem generally does not work well.

• Thank you again for taking the time. +1 – halirutan Jun 11 '18 at 8:44