High-order, Dispersionless "Fast-Hybrid" Wave Equation Solver. Part I: O(1)\mathcal{O}(1) Sampling Cost via Incident-Field Windowing and Recentering

This paper proposes a frequency/time hybrid integral-equation method for the time dependent wave equation in two and three-dimensional spatial domains. Relying on Fourier Transformation in time, the method utilizes a fixed (time-independent) number of frequency-domain integral-equation solutions to evaluate, with superalgebraically-small errors, time domain solutions for arbitrarily long times. The approach relies on two main elements, namely, 1) A smooth time-windowing methodology that enables accurate band-limited representations for arbitrarily-long time signals, and 2) A novel Fourier transform approach which, in a time-parallel manner and without causing spurious periodicity effects, delivers numerically dispersionless spectrally-accurate solutions. A similar hybrid technique can be obtained on the basis of Laplace transforms instead of Fourier transforms, but we do not consider the Laplace-based method in the present contribution. The algorithm can handle dispersive media, it can tackle complex physical structures, it enables parallelization in time in a straightforward manner, and it allows for time leaping---that is, solution sampling at any given time $T$ at $\mathcal{O}(1)$-bounded sampling cost, for arbitrarily large values of $T$, and without requirement of evaluation of the solution at intermediate times. The proposed frequency-time hybridization strategy, which generalizes to any linear partial differential equation in the time domain for which frequency-domain solutions can be obtained (including e.g. the time-domain Maxwell equations), and which is applicable in a wide range of scientific and engineering contexts, provides significant advantages over other available alternatives such as volumetric discretization, time-domain integral equations, and convolution-quadrature approaches.Comment: 33 pages, 8 figures, revised and extended manuscript (and now including direct comparisons to existing CQ and TDIE solver implementations) (Part I of II

Efficient evaluation of expectations of functions of a L\'evy process and its extremum

We prove simple general formulas for expectations of functions of a L\'evy process and its running extremum. Under additional conditions, we derive analytical formulas using the Fourier/Laplace inversion and Wiener-Hopf factorization, and discuss efficient numerical methods for realization of these formulas. As applications, the cumulative probability distribution function of the process and its running maximum and the price of the option to exchange the power of a stock for its maximum are calculated. The most efficient numerical methods use the sinh-acceleration technique and simplified trapezoid rule. The program in Matlab running on a Mac with moderate characteristics achieves the precision E-7 and better in several milliseconds, and E-14 - in a fraction of a second

A Wiener--Hopf Monte Carlo simulation technique for L\'{e}vy processes

We develop a completely new and straightforward method for simulating the joint law of the position and running maximum at a fixed time of a general L\'{e}vy process with a view to application in insurance and financial mathematics. Although different, our method takes lessons from Carr's so-called "Canadization" technique as well as Doney's method of stochastic bounds for L\'{e}vy processes; see Carr [Rev. Fin. Studies 11 (1998) 597--626] and Doney [Ann. Probab. 32 (2004) 1545-1552]. We rely fundamentally on the Wiener-Hopf decomposition for L\'{e}vy processes as well as taking advantage of recent developments in factorization techniques of the latter theory due to Vigon [Simplifiez vos L\'{e}vy en titillant la factorization de Wiener-Hopf (2002) Laboratoire de Math\'{e}matiques de L'INSA de Rouen] and Kuznetsov [Ann. Appl. Probab. 20 (2010) 1801--1830]. We illustrate our Wiener--Hopf Monte Carlo method on a number of different processes, including a new family of L\'{e}vy processes called hypergeometric L\'{e}vy processes. Moreover, we illustrate the robustness of working with a Wiener--Hopf decomposition with two extensions. The first extension shows that if one can successfully simulate for a given L\'{e}vy processes then one can successfully simulate for any independent sum of the latter process and a compound Poisson process. The second extension illustrates how one may produce a straightforward approximation for simulating the two-sided exit problem.Comment: Published in at http://dx.doi.org/10.1214/10-AAP746 the Annals of Applied Probability (http://www.imstat.org/aap/) by the Institute of Mathematical Statistics (http://www.imstat.org

On the evaluation of quasi-periodic Green functions and wave-scattering at and around Rayleigh-Wood anomalies

This article presents full-spectrum, well-conditioned, Green-function methodologies for evaluation of scattering by general periodic structures, which remains applicable on a set of challenging singular configurations, usually called Rayleigh-Wood (RW) anomalies (at which the quasi-periodic Green function ceases to exist), where most existing quasi-periodic solvers break down. After reviewing a variety of existing fast-converging numerical procedures commonly used to compute the classical quasi-periodic Green-function, the present work explores the difficulties they present around RW-anomalies and introduces the concept of hybrid “spatial/spectral” representations. Such expressions allow both the modification of existing methods to obtain convergence at RW-anomalies as well as the application of a slight generalization of the Woodbury-Sherman-Morrison formulae together with a limiting procedure to bypass the singularities. (Although, for definiteness, the overall approach is applied to the scalar (acoustic) wave-scattering problem in the frequency domain, the approach can be extended in a straightforward manner to the harmonic Maxwell's and elasticity equations.) Ultimately, this thorough study of RW-anomalies yields fast and highly-accurate solvers, which are demonstrated with a variety of simulations of wave-scattering phenomena by arrays of particles, crossed impenetrable and penetrable diffraction gratings and other related structures. In particular, the methods developed in this article can be used to “upgrade” classical approaches, resulting in algorithms that are applicable throughout the spectrum, and it provides new methods for cases where previous approaches are either costly or fail altogether. In particular, it is suggested that the proposed shifted Green function approach may provide the only viable alternative for treatment of three-dimensional high-frequency configurations with either one or two directions of periodicity. A variety of computational examples are presented which demonstrate the flexibility of the overall approach

Stable interpolation with exponential-polynomial splines and node selection via greedy algorithms

In this work we extend some ideas about greedy algorithms, which are well-established tools for, e.g., kernel bases, and exponential-polynomial splines whose main drawback consists in possible overfitting and consequent oscillations of the approximant. To partially overcome this issue, we develop some results on theoretically optimal interpolation points. Moreover, we introduce two algorithms which perform an adaptive selection of the spline interpolation points based on the minimization either of the sample residuals (f-greedy), or of an upper bound for the approximation error based on the spline Lebesgue function (λ-greedy). Both methods allow us to obtain an adaptive selection of the sampling points, i.e., the spline nodes. While the f-greedy selection is tailored to one specific target function, the λ-greedy algorithm enables us to define target-data-independent interpolation nodes