# Numerically computing the average value of an oscillating function

I want to compute the average value of a function $f(t)$, defined as $$\bar{f} = \lim_{T \rightarrow \infty} \frac{1}{T} \int_0^T f(t) \, \textrm{d}t \, .$$ where $f$ obtained numerically. The functions is cheap enough to not need to pre-compute and is bounded, but has quite fast oscillations. Numerically I have been approximating this by setting $T = T_0$ (large but finite) and evaluating the integral using scipy's QUADPACK methods, then increasing $T_0$ to test that the values $\bar{f}(T_0)$ converge.

It's made more difficult though because $f(t)$ is an oscillating function so fixing the end value of the integral changes the value you find. My first thought was to average over some range of $T$ around $T_0$ with a second integral $$\frac{1}{a} \int_{T_0}^{T_0+a} \textrm{d}T \left( \frac{1}{T} \int_0^T f(t) \, \textrm{d}t \right)$$ (where I would perhaps choose $a$ as 1% or 5% of $T_0$) but being a double integral this makes the QUADPACK evaluation much slower, and ignores the fact that the integrand is separable.

I realised after that a smarter version is to break the integration range into the rectangular part for $t \leq T_0$ and the triangular region for $t > T_0$. In the rectangular region the integrals over $t$ and $T$ are independent and I only need to compute the 1D integral over $t$ using QUADPACK. The $T$ integral is simple and can be done by hand. This reduces the area I need to perform the 2D integral over to just the triangular region.

Someone pointed out to me afterwards that this trick can be iterated by separating the triangle into a square (in which, by the same argument, we only need to perform a 1D integral) and two smaller triangles. This could be repeated again on the smaller triangles, but each time we do this the number of 2D integrals we need to perform doubles. Presumably after some number of iterations splitting the region further adds no benefit (and maybe even makes the computation more expensive), is this right?

Is there a way to estimate how the speed of QUADPACK's double integration scales with the size of the region so that I can argue when would be a good point to terminate this sub-division?

Are there better ways to approach the original problem?