Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Laplace equation in cylindrical coordinates

  1. Nov 6, 2009 #1
    Can anyone help with the solution of the Laplace equation in cylindrical coordinates
    [tex]\frac{\partial^{2} p}{\partial r^{2}}[/tex] [tex]+[/tex] [tex]\frac{1}{r}[/tex] [tex]\frac{\partial p}{\partial r}[/tex] [tex]+[/tex] [tex]\frac{\partial^{2} p}{\partial z^{2}}[/tex] [tex]= 0 [/tex]

    with Neumann no-flux boundaries:
    [tex]\frac{\partial p}{\partial r}[/tex] [tex]\left(0,z\right)[/tex] [tex] = 0[/tex]
    [tex]\frac{\partial p}{\partial r}[/tex] [tex]\left(Rmax,z\right)[/tex] [tex] = 0[/tex]
    [tex]\frac{\partial p}{\partial z}[/tex] [tex]\left(r,Zmax\right)[/tex] [tex] = 0[/tex]

    and a Dirichlet upper boundary:
    [tex]p(r,0) = f(r)[/tex]
     
  2. jcsd
  3. Nov 6, 2009 #2

    arildno

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member
    Dearly Missed

    With a trial solution using separation of variables, p(r,z)=P(r)*Q(z), we get:
    [tex]\frac{P''(r)+\frac{1}{r}P'(r)}{P(r)}=\frac{1}{Q(z)}Q''(z)[/tex]

    Thus, we get two ordinary diff.eqs,
    [tex]P''(r)+\frac{1}{r}P'(r)=CP(r) (*)[/tex]
    [tex]Q''(z)=CQ(z)(**)[/tex]

    You might try to work with these two, (*) clearly being the most problematic one.
     
  4. Nov 6, 2009 #3

    LCKurtz

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member

    But maybe not too problematic. Multiply it through by r2 and you have a form of Bessel's equation.
     
  5. Nov 7, 2009 #4
    Thanks for getting me started, and pardon the ignorance of this geochemist who hasn't taken a pde class. I'm trying to follow you, and can't understand how you got that 2nd ode. Shouldn't (**) be

    [tex]Q''(z) = -C [/tex] [tex]Q(z)[/tex] [tex] ?[/tex]
     
  6. Nov 7, 2009 #5

    arildno

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member
    Dearly Missed

    You are BOTH right. Mea culpa.
     
  7. Nov 7, 2009 #6
    Moving along here...

    So the solution to (**):

    [tex] Q(z) = Acosh(\lambda z)+Bsinh(\lambda z)[/tex]

    As Arildno predicted, the solution to (*) is more problematic, since this is my first experience with Bessel functions. Here's where I am so far...

    [tex]P(r) = CJ_{0}(\lambda r)+DY_{0}(\lambda r)[/tex]

    And here I set the second term must go to zero for the solution to be bounded.

    So the final solution is of the form:
    [tex]p(r,z) = [Acosh(\lambda z)+Bsinh(\lambda z)]J_{0}(\lambda r)[/tex]

    Using the 2nd Neumann BC:
    [tex]\frac{\partial p}{\partial r}(r=Rmax=9.5) = 0[/tex]
    This means that
    [tex]J'_{0}(9.5\lambda) = \frac{\partial J_{0}}{\partial r}= 0[/tex]
    But I'm not sure how this helps me right now, so I'll keep it in mind...

    Using the 3rd Neumann BC:
    [tex]\frac{\partial p}{\partial z}(z=Zmax=20) = 0[/tex]
    This means that
    [tex][Acosh(20\lambda)+Bsinh(20\lambda)]=0[/tex]
    Ah, this is more helpful. I learned that A=-B, and I think at the bottom boundary, this whole term drops out.

    Using the 1st Neumann BC: (I guess I could have done these in order...)
    [tex]\frac{\partial p}{\partial r}(r=0) = 0[/tex]
    This means that
    [tex]J'_{0}(0\lambda) = \frac{\partial J_{0}}{\partial r}= 0[/tex]
    Again I'm not really sure how this helps me right now, so I'll keep it in mind...

    So now we have:
    [tex]p(r,z)=A[cosh(\lambda z)-sinh(\lambda z)]J_{0}(\lambda r)[/tex]

    with
    [tex]J'_{0}(0\lambda, r=0) = J'_{0}(9.5\lambda, r=9.5) = 0[/tex]

    Now here's where it gets tricky and I get stuck...

    Using the fourth Dirichlet upper (z=0) boundary condition, I get:

    [tex]p(r,0)=f(r)=\sum_{\infty}^{n=1}A_{n}[cosh(0)-sinh(0)]J_{0}(\lambda _{n}r)=\sum_{\infty}^{n=1}A_{n}J_{0}(\lambda _{n}r)[/tex]

    Unfortunately, f(r) is not a constant, but rather an ugly empirical function of r.

    So I have written:

    [tex]f(r)=\sum_{\infty}^{n=1}A_{n}J_{0}(\lambda _{n}r)[/tex]

    And I have no idea where to go from here. I gather this will eventually end up in Matlab or R, but first I need a better understanding of the Bessel function J0 and A. I've read that the Bessel function of the first kind of order n (in my case n=0, right?) can be expressed as:

    [tex]J_{n}(x) = \sum_{\infty}^{k=1}\frac {(-1)^{k}(x/2)^{n+2k}}{k!\Gamma (n+k+1)}}[/tex]

    Again, any help would be greatly appreciated.
     
  8. Nov 7, 2009 #7

    arildno

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member
    Dearly Missed

    Now, assuming that the class of J_0's, with scaled arguments, represents a COMPLETE BASIS for functions on R (i.e, that any function is representable as a linear combination of these J_0's), then the A_n's are simply the required coefficients.
     
  9. Nov 7, 2009 #8

    LCKurtz

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member

    I haven't checked all your details, but your P equation and boundary conditions appear to be a Sturm-Liouville system which answers questions about orthogonality of the eigenfunctions, gives formulas for the eigenfunction expansion coefficients and settles convergence. For example, see:

    http://www.efunda.com/math/ode/Sturm_liouville.cfm
     
  10. Nov 7, 2009 #9
    Thanks again. Arildno, you mentioned that the [tex]A_{n}[/tex]'s are the required coefficients. But if I want to define a pressure distribution in this cylinder, don't I also need to figure out [tex]\lambda[/tex]? And is it true that [tex]\lambda[/tex] in P(r) and [tex]\lambda[/tex] in Q(z) are not the same?
     
  11. Nov 7, 2009 #10

    arildno

    User Avatar
    Science Advisor
    Homework Helper
    Gold Member
    Dearly Missed

    No, the [itex]\lambda_{n}[/itex] MUST be the same numbers; otherwise, your diff.eq won't be satisfied. (Remember, it is directly related to your C's!!)
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook