Evaluating matrix functions for exponential integrators via Carathéodory-Fejér approximation and contour integrals

Among the fastest methods for solving stiff PDE are exponential integrators, which require the evaluation of $f(A)$, where $A$ is a negative definite matrix and $f$ is the exponential function or one of the related ``$\varphi$ functions'' such as $\varphi_1(z) = (e^z-1)/z$. Building on previous work by Trefethen and Gutknecht, Gonchar and Rakhmanov, and Lu, we propose two methods for the fast evaluation of $f(A)$ that are especially useful when shifted systems $(A+zI)x=b$ can be solved efficiently, e.g. by a sparse direct solver. The first method method is based on best rational approximations to $f$ on the negative real axis computed via the Carathéodory-Fejér procedure, and we conjecture that the accuracy scales as $(9.28903\dots)^{-2n}$, where $n$ is the number of complex matrix solves. In particular, three matrix solves suffice to evaluate $f(A)$ to approximately six digits of accuracy. The second method is an application of the trapezoid rule on a Talbot-type contour

FIESTA 2: parallelizeable multiloop numerical calculations

The program FIESTA has been completely rewritten. Now it can be used not only as a tool to evaluate Feynman integrals numerically, but also to expand Feynman integrals automatically in limits of momenta and masses with the use of sector decompositions and Mellin-Barnes representations. Other important improvements to the code are complete parallelization (even to multiple computers), high-precision arithmetics (allowing to calculate integrals which were undoable before), new integrators and Speer sectors as a strategy, the possibility to evaluate more general parametric integrals.Comment: 31 pages, 5 figure

Efficient approximation of functions of some large matrices by partial fraction expansions

Some important applicative problems require the evaluation of functions $\Psi$ of large and sparse and/or \emph{localized} matrices $A$. Popular and interesting techniques for computing $\Psi(A)$ and $\Psi(A)\mathbf{v}$, where $\mathbf{v}$ is a vector, are based on partial fraction expansions. However, some of these techniques require solving several linear systems whose matrices differ from $A$ by a complex multiple of the identity matrix $I$ for computing $\Psi(A)\mathbf{v}$ or require inverting sequences of matrices with the same characteristics for computing $\Psi(A)$. Here we study the use and the convergence of a recent technique for generating sequences of incomplete factorizations of matrices in order to face with both these issues. The solution of the sequences of linear systems and approximate matrix inversions above can be computed efficiently provided that $A^{-1}$ shows certain decay properties. These strategies have good parallel potentialities. Our claims are confirmed by numerical tests

Algorithms to Evaluate Multiple Sums for Loop Computations

We present algorithms to evaluate two types of multiple sums, which appear in higher-order loop computations. We consider expansions of a generalized hypergeometric-type sums, \sum_{n_1,...,n_N} [Gamma(a1.n+c1) Gamma(a2.n}+c2) ... Gamma(aM.n+cM)] / [Gamma(b1.n+d1) Gamma(b2.n+d2) ... Gamma(bM.n+dM)] x1^n1...xN^nN with $ai.n=\sum_{j=1}^N a_{ij}nj$, etc., in a small parameter epsilon around rational values of ci,di's. Type I sum corresponds to the case where, in the limit epsilon -> 0, the summand reduces to a rational function of nj's times x1^n1...xN^nN; ci,di's can depend on an external integer index. Type II sum is a double sum (N=2), where ci,di's are half-integers or integers as epsilon -> 0 and xi=1; we consider some specific cases where at most six Gamma functions remain in the limit epsilon -> 0. The algorithms enable evaluations of arbitrary expansion coefficients in epsilon in terms of Z-sums and multiple polylogarithms (generalized multiple zeta values). We also present applications of these algorithms. In particular, Type I sums can be used to generate a new class of relations among generalized multiple zeta values. We provide a Mathematica package, in which these algorithms are implemented.Comment: 30 pages, 2 figures; address of Mathematica package in Sec.6; version to appear in J.Math.Phy

Lepton-pair scattering with an off-shell and an on-shell photon at two loops in massless QED

We compute the two-loop QED helicity amplitudes for the scattering of a lepton pair with an off-shell and an on-shell photon, $0\to\ell\bar\ell\gamma\gamma^*$, using the approximation of massless leptons. We express all master integrals relevant for the scattering of four massless particles with a single external off-shell leg up to two loops in a basis of algebraically independent multiple polylogarithms, which guarantees an efficient numerical evaluation and compact analytic representations of the amplitudes. Analytic forms of the amplitudes are reconstructed from numerical evaluations over finite fields. Our results complete the amplitude-level ingredients contributing to the N$^3$LO predictions of electron-muon scattering $e\mu\to e\mu$, which are required to meet the precision goal of the future MUonE experiment.Comment: 16 + 13 pages, 2 figures, 2 table

Lepton-pair scattering with an off-shell and an on-shell photon at two loops in massless QED

We compute the two-loop QED helicity amplitudes for the scattering of a lepton pair with an off-shell and an on-shell photon, 0 → ℓℓ¯γγ*, using the approximation of massless leptons. We express all master integrals relevant for the scattering of four massless particles with a single external off-shell leg up to two loops in a basis of algebraically independent multiple polylogarithm, which guarantees an efficient numerical evaluation and compact analytic representations of the amplitudes. Analytic forms of the amplitudes are reconstructed from numerical evaluations over finite fields. Our results complete the amplitude-level ingredients contributing to the N3LO predictions of electron-muon scattering eμ → eμ, which are required to meet the precision goal of the future MUonE experiment