[SciPy-User] [SciPy-user] Speeding up scipy.integrate.quad

[Warren Weckesser]

On Tue, Oct 25, 2011 at 4:48 AM, tomrichardson <
thomas.d.richardson at kcl.ac.uk> wrote:

>
> I'm currently trying to implement a Metropolis Hastings algorithm for which
> a
> numerical integration has to be called to obtain the acceptance probability
> of any given trial in the iteration. I'm currently using the
> scipy.integrate.quad function but it's too slow.Only part of the integrand
> has an analytical form - in fact it depends on an interpolation function
> that uses data in two numpy arrays. I've plotted the integrand and it
> diverges near the lower limit of integration, there is a term
>
> 1/(r^2-R^2) where r is the the dummy variable over which I am integrating
> and R is the lower limit of integration [I've already added a small
> increment to R to avoid a zero division error]
>


The integral of 1/(r^2 - R^2) from r=R to r=a>R is divergent (i.e. the area
under that curve is infinite).  Does some other term in your integrand
cancel the singularity at r=R to give a convergent integral?  If so, can you
rewrite your integrand so that these terms are handled separately (ideally
analytically), and only use quad for the well-behaved part of the integrand?

Warren



>
> Things I've tried:
> -I've implemented all of the other QUADPACK integration methods in scipy
> which are all slower
> -reduced the epsabs and epsrel limits
> -split the interval into subranges that weight the divergence
>
> Are there any special purpose integrators that anyone can recommend that
> can
> deal with this integrand and compete favourably with quad for speed?
>
> Ideas:
>
> I'm not familiar with the C language but I have Cython and Weave installed
> as part of the Enthought Python pack. Would implementing the numerical
> integration in the C compiler increase the speed dramatically? If so how do
> I introduce the numpy arrays in the cdef terminology for Cython / are there
> any C++ recipes that use a sampled data input that I can implement
> (relatively pain free as a C++ noob).
>
> --
> View this message in context:
> http://old.nabble.com/Speeding-up-scipy.integrate.quad-tp32716200p32716200.html
> Sent from the Scipy-User mailing list archive at Nabble.com.
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.scipy.org/pipermail/scipy-user/attachments/20111025/6b793255/attachment.html>