0
$\begingroup$

I need a fast and accurate method to calculate 3d spherical volume integrals.

I have pre-calculated data of high precision that just needs a few trivial manipulations before each integration step - in other words function calls are mostly just array look ups and multiplication and so are relatively cheap.

The data is evenly spaced in spherical polar coordinates: 4 samples per degree for theta and phi, 500 radial samples.

I need good precision, but I also need to perform several million such integrals, so speed is an important concern.

What algorithm is most suitable for this problem?

Bonus points for a link to an implementation in C++!

Improve this question
$\endgroup$
2
  • 1
    $\begingroup$ Do you need your data evenly spaced? Probably that's not the best way to perform the integrals. Have you checked on Lebedev quadrature? There is a C implementation in Github. $\endgroup$
    – nicoguaro ♦
    Jun 16, 2017 at 14:40
  • 1
    $\begingroup$ ... or, a bit simpler than the Lebedev stuff, use Legendre nodes for the polar angle. Equidistant in $\phi$ plus Legendre in $\theta$ are optimal by themselves, Lebedev adds the two-dimensional view. When I used it, I concluded for myself that the reduction in gridpoints is not worth the complication of using the Lebedev points). $\endgroup$
    – davidhigh
    Jun 16, 2017 at 15:01

2 Answers 2

Reset to default
2
$\begingroup$

So you want to evaluate

$$ F=\int_0^r \int_0^{2\pi} \int_{0}^\pi f(r,\phi,\theta) \, r^2 \, \sin (\theta) \ dr \, d\phi \, d\theta $$ and you have your function $f$ available on equidistant gridpoints $r_i$, $\phi_j$ and $\theta_k$. Call

$$ f_{ijk} = f(r_i, \phi_j, \theta_k) $$

The plain simple approach is then just to use the trapezoidal rule (I'll omit the one-halves at the endpoints):

$$ F \approx \sum_{ijk}f_{ijk} \, r_i^2 \, \sin(\theta_i) $$

... or, if you're feeling for higher order, use another Newton-Cotes formula.

Just to note that once. For other approaches, see the comment section.

EDIT: I've just noticed you were asking for the "best" method, so this will probably be quite boring for you.

Improve this answer
$\endgroup$
0
$\begingroup$

If the subvolumes are evenly distributed and you can assume the function is constant inside them you probably only need a lookup table for the volume size in the radius and do a weighted sum. If the assumptions are not true it will probably not be very accurate though :)

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.