[SciPy-user] Integrating Array

[rgold at lsw.uni-heidelberg.de]

First of all, thanks for the extensive and quick answers!
Unfortunately, I do not have an analytic expression of my integrand!
That's why I can't use quad! More precisely I perform the integral over a
product of two arrays from which one is the result of a function
calculating Legendre Polynomials, but the other array is pure non-analytic
data.
I have posted here a few months ago my problem when calculating the
Legendre-Polynomials using scipy.special.legendre() and thanks to your
help and the help of others I have written a function PLs(l,x) which
calculates for each given value x between -1 and 1 the value of the
Legendre Polynomial l and now everything is stable AND accurate (except
for my integration)!

Actually I DO KNOW more about my integrand (but only about the one part
with the Legendre Polynomials) thats why I started interpolating (@Anne:
You are completely right about the other part (=the data-array)! I just
have not enough points,end of story! I think that the error comes from the
Legendre Polys:-).
This morning I've seen that the interpolation indeed improves the accuracy
allthough I am using trapz(). Does anyone of you have an explanation for
the failure of the simps routine (happens only when I interpolate,works
fine for the non-interpolated data:-/)?
Today I am trying to interpolate such, that I can use romberg! Let's see
what will happen...

Thanks again!
Roman

> The reason the clever routines (e.g., scipy.integrate.quad) are able
> to produce high accuracy is because they are able to sample your
> function more densely in regions where it is more complicated. Between
> samples, they generally make some very simple assumption about the
> function's behaviour - a polynomial of low order, usually. If you have
> already sampled your function and are not willing to sample it
> further, there's not much that can be done to improve on simps. If
> your function were very smooth, you could try very high order
> polynomials, but this  is unlikely to help you much. The real problem
> is that you only know your function at some points, and interpolation
> is the best guess you have at what happens between them. The short
> answer is if your function oscillates between points, you need more
> points.
>
> The best and easiest choice is to simply use quad on your function,
> and let it decide where to sample. From your question, you presumably
> have a good reason not to do that. If it's that your evaluation
> routine is much more efficient when you give it a whole vector to
> evaluate at once, there are choices:
>
> * scipy.integrate.quadrature is based on Gaussian quadrature. The idea
> here is to evaluate your function at carefully-chosen abscissas, then
> use carefully-chosen weights to produce a quadrature rule. The
> abscissas and weights can be chosen so that you get exact integrals
> for polynomials up to a certain order times a weight function. Anyway,
> it's good, and it evaluates your function at a vector of points at
> once. scipy.integrate.quadrature evaluates it on denser and denser
> grids until the result appears to converge. (If you need special
> calling conventions for your function, you can look into the code of
> scipy.integrate.quadrature, it's not very complicated.)
>
> * If your function can only be evaluated at an evenly-spaced grid, you
> can use scipy.integrate.simps on finer and finer grids until it
> appears to have converged. If you have slightly more flexibility, you
> can subdivide individual grid steps separately (so that regions where
> the function is more complicated get evaluated more frequently).
>
> * On the remote chance that you are getting Fourier coefficients and
> wanting to integrate the time-domain function, some clever algebra
> involving the inverse Fourier transform lets you produce integrals
> directly (using an IFFT if desired).
>
> We may be more able to help you if you tell us more about your problem
> (for example why you're not using quad).
>
> Good luck,
> Anne
> _______________________________________________
> SciPy-user mailing list
> SciPy-user at scipy.org
> http://projects.scipy.org/mailman/listinfo/scipy-user
>