A practical algorithm for computing Cauchy principal value integrals of oscillatory functions (Q426517)

The paper describes a numerical method for the approximation of Cauchy principal value integrals with oscillatory integrands, \[ \int_{-1}^1 (x-\tau)^{-1} f(x) \exp(i \omega x) dx, \] with some \(\tau \in (-1,1)\). The algorithm proceeds as follows. First, approximate the function \(f\) in the integrand by an interpolating polynomial \(p_N\) of degree \(N\), using the Clenshaw-Curtis nodes \(x_j = -\cos j \pi /N\), \(j = 0,1,\dots, N\). The polynomial \(p_N\) is expressed in terms of Chebyshev polynomials. Then, a simple and cheap approximation \(R\) is computed such that \(R \approx p_N(\tau)\), and the approximate integral \(\int_{-1}^1 (x-\tau)^{-1} p_N(x) \exp(i \omega x) dx\) is rewritten by adding and subtracting \(R\) in the numerator. In this respect, the approach is closely related to the ideas proposed by the reviewer [``Gaussian quadrature formulae of the third kind for Cauchy principal value integrals: Basic properties and error estimates'', J. Comput.\ Appl.\ Math.\ 65, No.\ 1--3, 97--114 (1995; Zbl 0847.41019)]. The indefinite Cauchy principal value integral \[ \int (p_N(x) - R) \exp(i \omega x) (x-\tau)^{-1} dx \] (that would be a regular integral if \(R = p_N(\tau)\)) is then also approximately expanded in the form of a product of \(\exp(i \omega x)\) and a sum of Chebyshev polynomials of degree up to \(M = N(1 + o(1))\). It is then possible to find a computationally cheap way to compute the coefficients of this expansion. A combination with the original rewriting of the approximate integral then yields the final quadrature formula. The author provides a detailed derivation of the method and a thorough error analysis, also taking into consideration the behaviour of the method as \(\tau\) and \(\omega\) are varied.