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

A 2D Finite Difference formulation in polar coordinates.

  1. May 2, 2017 #1
    So I have this PDE:

    d2T/dr2 + 1/r dT/dr + d2T/dθ2 = 0.

    How do I implement dT/dr || [r = 0] = 0? Also, what should I do about 1/r?

    This is actually the first time I am going to attack FDF in polar/cylindrical coordinates. I can finite-difference the base equation fairly decently; I am just having a hard time in implementing the derivative boundary condition at r = 0. Can someone give me an idea what to do?
  2. jcsd
  3. May 2, 2017 #2
    In practice, no computations are made on the line r=0. I assume you know how to discretize and how to obtain your finite difference matrix-vector system.
    Suppose you discretize and you have N+1 nodes from j=0..N in the radial direction where the nodes with index j=0 correspond to r=0. This leads to a matrix-vector system of size (N+1). The Neumann boundary condition (first order) can be written as ## \frac{T_1 - T_0}{\Delta R_1} = 0## or: ## T_0 = T_1##

    The solution on the computational nodes at r=0 are known when ##T_1## is known. So you replace every occurence of ##T_0## with ##T_1## and you delete the first row and column of the discretization matrix and you solve the size N matrix system. When the solution ## T_1 .. T_N## is found, you simply add ## T_0 = T_1## to your solution vector.
  4. May 3, 2017 #3
    Hi. Thanks for replying.

    I managed to do what you suggested in a 1-D system (in radial directions only). I wanted to stencil my problem in a 2-D system (radial and angular directions). This is my problem:


    I know the number of points is... :DD:DD:DD I'm doing this to kill off certain doubts.

    I stencil around T0, by using an energy balance approach. I got something like T0 = (2T1 + 2T3) / 4. Is this correct? This is actually the problem I am facing.

    If I implement the boundary condition that dT/dr [r = 0] = 0, that would mean either (T1-T0)/2Δr or (T3-T0)/2Δr is 0 right? Which do I replace? T1 or T3?
    Last edited: May 3, 2017
  5. May 3, 2017 #4


    User Avatar
    Homework Helper

    Physically [itex]r = 0[/itex] is an interior point of the domain. You can't impose a condition on the value of [itex]T[/itex] or [itex]\partial T/\partial r[/itex] there.
  6. May 3, 2017 #5
    If he is simulating only a quarter of a circle, then r=0 is on the boundary of the domain (it is not so clear in the text, but the problem description says ##\Delta \phi = \pi / 2## ).

    Then to implement the Neumann boundary condition:
    on the horizontal line, T1=T0
    on the 45-degree line, T2=T0
    on the vertical line, T3=T0
    So T0=T1=T2=T3
Share this great discussion with others via Reddit, Google+, Twitter, or Facebook

Have something to add?
Draft saved Draft deleted