Libsharp - spherical harmonic transforms revisited

We present libsharp, a code library for spherical harmonic transforms (SHTs), which evolved from the libpsht library, addressing several of its shortcomings, such as adding MPI support for distributed memory systems and SHTs of fields with arbitrary spin, but also supporting new developments in CPU instruction sets like the Advanced Vector Extensions (AVX) or fused multiply-accumulate (FMA) instructions. The library is implemented in portable C99 and provides an interface that can be easily accessed from other programming languages such as C++, Fortran, Python etc. Generally, libsharp's performance is at least on par with that of its predecessor; however, significant improvements were made to the algorithms for scalar SHTs, which are roughly twice as fast when using the same CPU capabilities. The library is available at http://sourceforge.net/projects/libsharp/ under the terms of the GNU General Public License

On the computation of Gaussian quadrature rules for Chebyshev sets of linearly independent functions

We consider the computation of quadrature rules that are exact for a Chebyshev set of linearly independent functions on an interval $[a,b]$. A general theory of Chebyshev sets guarantees the existence of rules with a Gaussian property, in the sense that $2l$ basis functions can be integrated exactly with just $l$ points and weights. Moreover, all weights are positive and the points lie inside the interval $[a,b]$. However, the points are not the roots of an orthogonal polynomial or any other known special function as in the case of regular Gaussian quadrature. The rules are characterized by a nonlinear system of equations, and earlier numerical methods have mostly focused on finding suitable starting values for a Newton iteration to solve this system. In this paper we describe an alternative scheme that is robust and generally applicable for so-called complete Chebyshev sets. These are ordered Chebyshev sets where the first $k$ elements also form a Chebyshev set for each $k$. The points of the quadrature rule are computed one by one, increasing exactness of the rule in each step. Each step reduces to finding the unique root of a univariate and monotonic function. As such, the scheme of this paper is guaranteed to succeed. The quadrature rules are of interest for integrals with non-smooth integrands that are not well approximated by polynomials

Quadrature Strategies for Constructing Polynomial Approximations

Finding suitable points for multivariate polynomial interpolation and approximation is a challenging task. Yet, despite this challenge, there has been tremendous research dedicated to this singular cause. In this paper, we begin by reviewing classical methods for finding suitable quadrature points for polynomial approximation in both the univariate and multivariate setting. Then, we categorize recent advances into those that propose a new sampling approach and those centered on an optimization strategy. The sampling approaches yield a favorable discretization of the domain, while the optimization methods pick a subset of the discretized samples that minimize certain objectives. While not all strategies follow this two-stage approach, most do. Sampling techniques covered include subsampling quadratures, Christoffel, induced and Monte Carlo methods. Optimization methods discussed range from linear programming ideas and Newton's method to greedy procedures from numerical linear algebra. Our exposition is aided by examples that implement some of the aforementioned strategies

Caratheodory-Tchakaloff Subsampling

We present a brief survey on the compression of discrete measures by Caratheodory-Tchakaloff Subsampling, its implementation by Linear or Quadratic Programming and the application to multivariate polynomial Least Squares. We also give an algorithm that computes the corresponding Caratheodory-Tchakaloff (CATCH) points and weights for polynomial spaces on compact sets and manifolds in 2D and 3D

Discontinuous collocation methods and gravitational self-force applications

Numerical simulations of extereme mass ratio inspirals, the mostimportant sources for the LISA detector, face several computational challenges. We present a new approach to evolving partial differential equations occurring in black hole perturbation theory and calculations of the self-force acting on point particles orbiting supermassive black holes. Such equations are distributionally sourced, and standard numerical methods, such as finite-difference or spectral methods, face difficulties associated with approximating discontinuous functions. However, in the self-force problem we typically have access to full a-priori information about the local structure of the discontinuity at the particle. Using this information, we show that high-order accuracy can be recovered by adding to the Lagrange interpolation formula a linear combination of certain jump amplitudes. We construct discontinuous spatial and temporal discretizations by operating on the corrected Lagrange formula. In a method-of-lines framework, this provides a simple and efficient method of solving time-dependent partial differential equations, without loss of accuracy near moving singularities or discontinuities. This method is well-suited for the problem of time-domain reconstruction of the metric perturbation via the Teukolsky or Regge-Wheeler-Zerilli formalisms. Parallel implementations on modern CPU and GPU architectures are discussed.Comment: 29 pages, 5 figure

WavePacket: A Matlab package for numerical quantum dynamics. I: Closed quantum systems and discrete variable representations

WavePacket is an open-source program package for the numerical simulation of quantum-mechanical dynamics. It can be used to solve time-independent or time-dependent linear Schr\"odinger and Liouville-von Neumann-equations in one or more dimensions. Also coupled equations can be treated, which allows to simulate molecular quantum dynamics beyond the Born-Oppenheimer approximation. Optionally accounting for the interaction with external electric fields within the semiclassical dipole approximation, WavePacket can be used to simulate experiments involving tailored light pulses in photo-induced physics or chemistry.The graphical capabilities allow visualization of quantum dynamics 'on the fly', including Wigner phase space representations. Being easy to use and highly versatile, WavePacket is well suited for the teaching of quantum mechanics as well as for research projects in atomic, molecular and optical physics or in physical or theoretical chemistry.The present Part I deals with the description of closed quantum systems in terms of Schr\"odinger equations. The emphasis is on discrete variable representations for spatial discretization as well as various techniques for temporal discretization.The upcoming Part II will focus on open quantum systems and dimension reduction; it also describes the codes for optimal control of quantum dynamics.The present work introduces the MATLAB version of WavePacket 5.2.1 which is hosted at the Sourceforge platform, where extensive Wiki-documentation as well as worked-out demonstration examples can be found