Driscoll&Healy integration for Multipole

Create issue
Issue #215 closed
Erik Schnetter created an issue

The attached patch implements a more accurate integration over the sphere, using an algorithm by Driscoll & Healy. This algorithm uses Gaussian integration weights, leading (almost) to exponential convergence.

Keyword: Multipole

Comments (8)

  1. Roland Haas
    • removed comment

    this seems to be buggy. The results differ by what looks like Pi from the results of test_simpson.par which is otherwise an identical parfile. Testing this with a constant value fr[i]=1.0*sin(th) in utils.cc also yields an answer 3.14 times too large when driscoll&Healy is used as compared to Simpson (which gives the correct 4*Pi result).

    I don't know how to fix this since I have no idea where the extra Pi creeps in.

    Do not apply in this form. :-)

  2. Erik Schnetter reporter
    • changed status to resolved
    • removed comment

    I found the missing factor pi. This factor is actually present in the (unnumbered!) equations in the Driscoll&Healy paper. In my tests, I left it out because I implicitly assumed that theta ranges [0...1], whereas of course it ranger [0...pi]. This made the results too large by a factor of pi.

    There was also an error in the index bounds. Driscoll&Healy want a point on the equator, which corresponds to even ntheta, not odd ntheta as previously in the code.

    Applied.

  3. Log in to comment