Quadrature methods for highly oscillatory singular integrals

We address the evaluation of highly oscillatory integrals, with power-law and logarithmic singularities. Such problems arise in numerical methods in engineering. Notably, the evaluation of oscillatory integrals dominates the run-time for wave-enriched boundary integral formulations for wave scattering, and many of these exhibit singularities. We show that the asymptotic behaviour of the integral depends on the integrand and its derivatives at the singular point of the integrand, the stationary points and the endpoints of the integral. A truncated asymptotic expansion achieves an error that decays faster for increasing frequency. Based on the asymptotic analysis, a Filon-type method is constructed to approximate the integral. Unlike an asymptotic expansion, the Filon method achieves high accuracy for both small and large frequency. Complex-valued quadrature involves interpolation at the zeros of polynomials orthogonal to a complex weight function. Numerical results indicate that the complex-valued Gaussian quadrature achieves the highest accuracy when the three methods are compared. However, while it achieves higher accuracy for the same number of function evaluations, it requires significant additional cost of computation of orthogonal polynomials and their zeros

Asymptotic expansions and fast computation of oscillatory Hilbert transforms

In this paper, we study the asymptotics and fast computation of the one-sided oscillatory Hilbert transforms of the form $$H^{+}(f(t)e^{i\omega t})(x)=-int_{0}^{\infty}e^{i\omega t}\frac{f(t)}{t-x}dt,\qquad \omega>0,\qquad x\geq 0,$$ where the bar indicates the Cauchy principal value and $f$ is a real-valued function with analytic continuation in the first quadrant, except possibly a branch point of algebraic type at the origin. When $x=0$, the integral is interpreted as a Hadamard finite-part integral, provided it is divergent. Asymptotic expansions in inverse powers of $\omega$ are derived for each fixed $x\geq 0$, which clarify the large $\omega$ behavior of this transform. We then present efficient and affordable approaches for numerical evaluation of such oscillatory transforms. Depending on the position of $x$, we classify our discussion into three regimes, namely, $x=\mathcal{O}(1)$ or $x\gg1$, $0<x\ll 1$ and $x=0$. Numerical experiments show that the convergence of the proposed methods greatly improve when the frequency $\omega$ increases. Some extensions to oscillatory Hilbert transforms with Bessel oscillators are briefly discussed as well.Comment: 32 pages, 6 figures, 4 table

Fast Computation of Fourier Integral Operators

We introduce a general purpose algorithm for rapidly computing certain types of oscillatory integrals which frequently arise in problems connected to wave propagation and general hyperbolic equations. The problem is to evaluate numerically a so-called Fourier integral operator (FIO) of the form $\int e^{2\pi i \Phi(x,\xi)} a(x,\xi) \hat{f}(\xi) \mathrm{d}\xi$ at points given on a Cartesian grid. Here, $\xi$ is a frequency variable, $\hat f(\xi)$ is the Fourier transform of the input $f$, $a(x,\xi)$ is an amplitude and $\Phi(x,\xi)$ is a phase function, which is typically as large as $|\xi|$; hence the integral is highly oscillatory at high frequencies. Because an FIO is a dense matrix, a naive matrix vector product with an input given on a Cartesian grid of size $N$ by $N$ would require $O(N^4)$ operations. This paper develops a new numerical algorithm which requires $O(N^{2.5} \log N)$ operations, and as low as $O(\sqrt{N})$ in storage space. It operates by localizing the integral over polar wedges with small angular aperture in the frequency plane. On each wedge, the algorithm factorizes the kernel $e^{2 \pi i \Phi(x,\xi)} a(x,\xi)$ into two components: 1) a diffeomorphism which is handled by means of a nonuniform FFT and 2) a residual factor which is handled by numerical separation of the spatial and frequency variables. The key to the complexity and accuracy estimates is that the separation rank of the residual kernel is \emph{provably independent of the problem size}. Several numerical examples demonstrate the efficiency and accuracy of the proposed methodology. We also discuss the potential of our ideas for various applications such as reflection seismology.Comment: 31 pages, 3 figure

Windowed Green Function Method for Nonuniform Open-Waveguide Problems

This contribution presents a novel Windowed Green Function (WGF) method for the solution of problems of wave propagation, scattering and radiation for structures which include open (dielectric) waveguides, waveguide junctions, as well as launching and/or termination sites and other nonuniformities. Based on use of a "slow-rise" smooth-windowing technique in conjunction with free-space Green functions and associated integral representations, the proposed approach produces numerical solutions with errors that decrease faster than any negative power of the window size. The proposed methodology bypasses some of the most significant challenges associated with waveguide simulation. In particular the WGF approach handles spatially-infinite dielectric waveguide structures without recourse to absorbing boundary conditions, it facilitates proper treatment of complex geometries, and it seamlessly incorporates the open-waveguide character and associated radiation conditions inherent in the problem under consideration. The overall WGF approach is demonstrated in this paper by means of a variety of numerical results for two-dimensional open-waveguide termination, launching and junction problems.Comment: 16 Page

Boundary integral methods in high frequency scattering

In this article we review recent progress on the design, analysis and implementation of numerical-asymptotic boundary integral methods for the computation of frequency-domain acoustic scattering in a homogeneous unbounded medium by a bounded obstacle. The main aim of the methods is to allow computation of scattering at arbitrarily high frequency with finite computational resources