# Evaluating oscillatory integrals with many independent periods and no closed forms

Most methods for oscillatory integrals I know about deal with integrals of the form $$\int f(x)e^{i\omega x}\,dx$$ where $\omega$ is large.

If I have an integral of the form $$\int f(x)g_1(x)\cdots g_n(x)\,dx,$$ where $g_k$ are oscillatory functions whose roots are only known approximately, but some kind of asymptotic form $$g_k(x) \sim e^{i\omega_k x}$$ is known, with the frequencies $\omega_k$ all different (and $\mathbb{Q}$-linearly independent), then how can I evaluate this integral?

Unlike in the case of $e^{i\omega x}$, the polynomial integrals $\int x^a \prod g_k(x)$ are not known, so I can't construct a set of polynomial interpolants for $f(x)$ and integrate the interpolants exactly.

In my exact problem, $g_k$'s are Bessel functions $J_0(\omega_k x)$, and $f(x)=x^\alpha$, and the region of integration is $[0,\infty)$. The method I am using now is to sum up integral contributions over intervals $[x_{k-1},x_k]$ between roots up to some cutoff $M$, then use the asymptotic expansion for $g_k(x)$ for large $x$. This algorithm's time complexity is exponential in $n$ because it involves expanding the product $g_1\ldots g_n$, each of which has a number $r$ of asymptotic terms, giving $r^n$ total terms; pruning terms that are too small doesn't reduce the run time enough to make this feasible for large $n$.

Heuristic non-rigorous answers, suggestions and references are all welcome.

I have worked on simpler integrals where there are points of stationary phase. I found two methods that work quite well.

One is to introduce an exponential damping factor that depends on the phase function, a kind of artificial viscosity if you like.

Another technique (where there are multiple points of stat. phase) was described in:

Tuck, E.O, Collins, J.L. and Wells, W.H., "On ship waves and their spectra", Journal of Ship Research, pp. 11–21, 1971.

That method applies exponential decay factors to the integrand where it becomes rapidly-oscillating away from the stat. phase points, but leaves the integrand intact where it is not.

That's me out of ideas!

• Thank you, but I don't quite see how this would work in this case. For one thing, there are no points of stationary phase on the real line, and the contributions from oscillations are significant to the final value, so must not be damped. – Kirill Jun 12 '15 at 22:14

As long as you have accurate values for the roots (or extrema) of the oscillatory part of your integrand, Longman's method (as I described in this answer) remains applicable. All you have to do is evaluate a bunch of integrals with intervals in between the roots using your favorite quadrature method, and treat these integrals as the terms of some alternating series. You can then use any number of convergence acceleration methods (Euler, Levin, Weniger, etc.) to "sum" this alternating series.

As an example, in this math.SE answer, I evaluated an infinite integral whose oscillatory part is a product of two Bessel functions.

• Wouldn't it matter that the roots are irregularly spaced (all of the periods are irrational and independent)? Why would you trust convergence acceleration for such an irregular sequence? – Kirill Feb 23 '16 at 12:20
• This was a while ago, I wanted to evaluate the integral to a thousand digits and if I remember correctly oscillatory quadrature was actually the first thing I tried. I don't remember the results, but I don't think it worked well at the time. – Kirill Feb 23 '16 at 12:28
• "Why would you trust convergence acceleration for such an irregular sequence?" - I wouldn't trust just one accelerator, tho. But, if at least three different accelerators are giving me consistent results, I'd think that the digits I got are at least plausible. FWIW, I've used Longman for infinite integrals of products of Bessel functions, and I've never been disappointed, especially when using Weniger's transformation as the accelerator. – J. M. Feb 23 '16 at 12:33
• The method I describe in the question is also an oscillatory quadrature method: expand the integrand in a series of terms of the form $x^a e^{b x}$, the infinite integrals for which have a closed form. I would trust such a method more than convergence acceleration. My understanding is they require something like strong monotonicity or a good understanding of error terms to be sure to work well. – Kirill Feb 23 '16 at 12:45
• If you can do a (generalized) Fourier expansion, then sure. – J. M. Feb 23 '16 at 12:51