0

I have a function and I will need to estimate its definite integral over many different intervals, some of them with endpoints quite close to each other,

\int_a^{b_k} f(x) dx

I could call integrate.quad afresh for each interval, but that would seem to be quite inefficicent as integrate.quad will be evaluating the function many times for each call. Is there any way to reuse the points that quad has already evaluated to make subsequent calls to quad faster?

(If not with integrate.quad, perhaps another scipy function or using a different public python library?)

Or perhaps I should cache my previous integral values and just call quad to compute the difference between the new point and the closest previously computed point?

Note also that the chosen endpoints are not all known in advance, they are chosen adaptively.

Nico Schlömer
  • 53,797
  • 27
  • 201
  • 249
Galen
  • 13
  • 3
  • `quad` can provide a lot of extra information (which you can read about just as well as I can), but I don't see a way of returning (or providing) actual integration points and values. But how expensive is it to evaluate your function? Caching and modifying values for the new case might be just as expensive as a 'raw' calculation. – hpaulj Jun 26 '21 at 18:17
  • 2
    Maybe you could use a different integration method. Supposedly an advantage of `quad` is it chooses an optimal (or at least good) number of evaluation points (I haven't studied the details). On the other hand with something like `trapezoid` you provide the array of `x` and `y` values. Depending on the function and required accuracy that may require more points than `quad`, but they can be evaluated "all-at-once" with `numpy` array methods. – hpaulj Jun 26 '21 at 18:22
  • Thanks. Evaluation is potentially rather expensive: O(1000) trig function evaluations, and I need it in real time. I found `quad` can return the subintervals and integral approximations on those subintervals. I'm not sure what method is used to compute the values on the subintervals (I verified it is not just trapezoid) but they do sum to the value of the complete integral. Maybe I can do the complete integral once, keep those points and their values, and then do my own interpolation to points between them. – Galen Jun 27 '21 at 01:14
  • The sequence `bₖ`, is it monotonic? – gboffi Jul 09 '21 at 09:43
  • Not necessarily, but a partial solution assuming monotonic b_k would still be useful. – Galen Jul 09 '21 at 19:15

1 Answers1

0

There is little hope that you can use general adaptive integration methods for intervals with different end points. The reason for this is that, in general, an integral value can change greatly if you change the limits by just a bit. You would have to somehow provide information that your function isn't "crazy" at the boundary, give estimates etc., and I know of no such method.

What you can do however is transform your integrals to the same interval, say [0, 1], by variable substitution. You end up with the task of integrating many (parametrized) functions over [0, 1].

\int_0^1 f_a(x)

You can then try and vectorize the function calls, i.e., evaluate all functions at once for a point or a set of points.

quadpy (a project of mine) supports vectorized quadrature which can dramatically speed up the integration. Example with x ** a for a number of different exponents:

import numpy as np
import quadpy

a = np.arange(1, 7)

def f(x):
    return np.power.outer(x, a).T


val, err = quadpy.quad(f, 0.0, 1.0)

print(val)
print(err)

[0.5        0.33333333 0.25       0.2        0.16666667 0.14285714]
[1.87519478e-20 2.58599737e-20 4.45666434e-20 5.19063171e-20
 5.45158477e-20 4.09060649e-20]

Note that in f, instead of iterating over the a the computation is vectorized. This greatly speeds up the computation.

Nico Schlömer
  • 53,797
  • 27
  • 201
  • 249
  • It looks like this is similar to what I need, but evaluating for many values of a parameter at once, not evaluating at many different interval end points. Also unfortunately my use case is adaptive, meaning the interval end point of subsequent evaluations may depend on the integreal at previously computed points, so I wouldn't be able to just vectorize the evaluation for all end points in any case. – Galen Jul 07 '21 at 17:14
  • Perhaps you should specify more precisely in your question what your problem is. – Nico Schlömer Jul 07 '21 at 17:54
  • You mean what exactly is the function I am integrating? I will need to do this with several functions that vary depending on the application, so I can't say specifically. I hope it is clear that my general problem is being able to efficiently evaluate multiple definite integrals of a single function where the endpoints may vary. – Galen Jul 08 '21 at 18:30
  • No, I mean your description of the problem. You write"I have a function and I will need to estimate its integral at many different points". What does this mean, "estimate at different points"? – Nico Schlömer Jul 08 '21 at 18:36
  • 1
    Thank you, I edited the original description to clarify: I have a function and I will need to estimate its definite integral over many different intervals, some of them with endpoints quite close to each other. – Galen Jul 09 '21 at 00:41
  • @Galen Edited my answer, perhaps it's useful to you. – Nico Schlömer Jul 09 '21 at 09:37