Search

Numerical integration over the real line for analytic functions is studied. Our main focus is on the sharpness of the error bounds. We first derive two general lower estimates for the worst-case integration error, and then apply these to establish lower bounds for various quadrature rules. These bounds turn out to either be novel or improve upon existing results, leading to lower bounds that closely match upper bounds for various formulas. Specifically, for the suitably truncated trapezoidal rule, we improve upon general lower bounds on the worst-case error obtained by Sugihara [Numer. Math., 75 (1997), pp. 379-395] and provide exceptionally sharp lower bounds apart from a polynomial factor, and in particular we show that the worst-case error for the trapezoidal rule by Sugihara is not improvable by more than a polynomial factor. Additionally, our research reveals a discrepancy between the error decay of the trapezoidal rule and Sugihara's lower bound for general numerical integration rules, introducing a new open problem. Moreover, the Gauss-Hermite quadrature is proven suboptimal under the decay conditions on integrands we consider, a result not deducible from upper-bound arguments alone. Furthermore, to establish the near-optimality of the suitably scaled Gauss-Legendre and Clenshaw-Curtis quadratures, we generalize a recent result of Trefethen [SIAM Rev., 64 (2022), pp. 132-150] for the upper error bounds in terms of the decay conditions. [ABSTRACT FROM AUTHOR]

We investigate a class of parametric elliptic eigenvalue problems with homogeneous essential boundary conditions, where the coefficients (and hence the solution) may depend on a parameter. For the efficient approximate evaluation of parameter sensitivities of the first eigenpairs on the entire parameter space we propose and analyze Gevrey class and analytic regularity of the solution with respect to the parameters. This is made possible by a novel proof technique, which we introduce and demonstrate in this paper. Our regularity result has immediate implications for convergence of various numerical schemes for parametric elliptic eigenvalue problems, in particular, for elliptic eigenvalue problems with infinitely many parameters arising from elliptic differential operators with random coefficients, e.g., integration by quasi--Monte Carlo methods. [ABSTRACT FROM AUTHOR]

We introduce a predictor-corrector discretization scheme for the numerical integration of a class of stochastic differential equations and prove that it converges with weak order 1.0. The key feature of the new scheme is that it builds up sequentially (and recursively) in the dimension of the state space of the solution, hence making it suitable for approximations of high-dimensional state space models. We show, using the stochastic Lorenz 96 system as a test model, that the proposed method can operate with larger time steps than the standard Euler--Maruyama scheme and, therefore, generate valid approximations with a smaller computational cost. We also introduce the theoretical analysis of the error incurred by the new predictor-corrector scheme when used as a building block for discrete-time Bayesian filters for continuous-time systems. Finally, we assess the performance of several ensemble Kalman filters that incorporate the proposed sequential predictor-corrector Euler scheme and the standard Euler--Maruyama method. The numerical experiments show that the filters employing the new sequential scheme can operate with larger time steps, smaller Monte Carlo ensembles, and noisier systems. [ABSTRACT FROM AUTHOR]

Partial differential equations (PDEs) with uncertain or random inputs have been considered in many studies of uncertainty quantification. In forward uncertainty quantification, one is interested in analyzing the stochastic response of the PDE subject to input uncertainty, which usually involves solving high-dimensional integrals of the PDE output over a sequence of stochastic variables. In practical computations, one typically needs to discretize the problem in several ways: approximating an infinite-dimensional input random field with a finite-dimensional random field, spatial discretization of the PDE using, e.g., finite elements, and approximating high-dimensional integrals using cubatures such as quasi--Monte Carlo methods. In this paper, we focus on the error resulting from dimension truncation of an input random field. We show how Taylor series can be used to derive theoretical dimension truncation rates for a wide class of problems and we provide a simple checklist of conditions that a parametric mathematical model needs to satisfy in order for our dimension truncation error bound to hold. Some of the novel features of our approach include that our results are applicable to nonaffine parametric operator equations, dimensionally truncated conforming finite element discretized solutions of parametric PDEs, and even compositions of PDE solutions with smooth nonlinear quantities of interest. As a specific application of our method, we derive an improved dimension truncation error bound for elliptic PDEs with lognormally parameterized diffusion coefficients. Numerical examples support our theoretical findings. [ABSTRACT FROM AUTHOR]

We introduce a numerical approach to computing the SchrÖdinger map (SM) based on the Hasimoto transform which relates the SM flow to a cubic nonlinear Schr\"odinger (NLS) equation. In exploiting this nonlinear transform we are able to introduce the first fully explicit unconditionally stable symmetric integrators for the SM equation. Our approach consists of two parts: an integration of the NLS equation followed by the numerical evaluation of the Hasimoto transform. Motivated by the desire to study rough solutions to the SM equation, we also introduce a new symmetric low-regularity integrator for the NLS equation. This is combined with our novel fast low-regularity Hasimoto (FLowRH) transform, based on a tailored analysis of the resonance structures in the Magnus expansion and a fast realization based on block-Toeplitz partitions, to yield an efficient low-regularity integrator for the SM equation. This scheme in particular allows us to obtain approximations to the SM in a more general regime (i.e., under lower-regularity assumptions) than previously proposed methods. The favorable properties of our methods are exhibited both in theoretical convergence analysis and in numerical experiments. [ABSTRACT FROM AUTHOR]

We consider a Markov chain Monte Carlo method, known as Hamiltonian Monte Carlo (HMC), on compact manifolds in Euclidean space. It utilizes Hamiltonian dynamics to generate samples approximating a target distribution in high dimensions efficiently. The efficiency of HMC is characterized by its convergence property, called geometric ergodicity. This property is important to generate low-correlated samples. It also plays a crucial role in establishing the error estimate for the quadrature of bounded functions by HMC sampling, referred to as the Hoeffding-type inequality. While the geometric ergodicity has been proved for HMC on Euclidean space, it has not been established on manifolds. In this paper, we prove the geometric ergodicity for HMC on compact manifolds. As an example to confirm the efficiency of the proposed HMC method, we consider a sampling problem associated with the N-vortex problem on the unit sphere, which is a statistical model of two-dimensional turbulence. We apply HMC to approximate the statistical quantities with respect to the invariant measure of the N-vortex problem, called the Gibbs measure. We observe the organization of large vortex structures as seen in two-dimensional turbulence. [ABSTRACT FROM AUTHOR]

A basic leapfrog integrator and its energy-preserving and variational/symplectic variants are proposed and studied for the numerical integration of the equations of motion of relativistic charged particles in an electromagnetic field. The methods are based on a four-dimensional formulation of the equations of motion. Structure-preserving properties of the numerical methods are analyzed, in particular conservation and long-time near-conservation of energy and mass shell as well as preservation of volume in phase space. In the non-relativistic limit, the considered methods reduce to the Boris algorithm for non-relativistic charged-particle dynamics and its energy-preserving and variational/symplectic variants. [ABSTRACT FROM AUTHOR]

A framework is presented to design multirate time stepping algorithms for two dissipative models with coupling across a physical interface. The coupling takes the form of boundary conditions imposed on the interface, relating the solution variables for both models to each other. The multirate aspect arises when numerical time integration is performed with different time step sizes for the component models. In this paper, we seek to identify a unified approach to develop multirate algorithms for these coupled problems. This effort is pursued though the use of discontinuous-Galerkin time stepping methods, acting as a general unified framework, with different time step sizes. The subproblems are coupled across user-defined intervals of time, called coupling windows, using polynomials that are continuous on the window. The coupling method is shown to reproduce the correct interfacial energy dissipation, discrete conservation of fluxes, and asymptotic accuracy. In principle, methods of arbitrary order are possible. As a first step, herein we focus on the presentation and analysis of monolithic methods for advection-diffusion models coupled via generalized Robintype conditions. The monolithic methods could be computed using a Schur-complement approach. The framework and analysis herein accommodate multirate coupling for a wide variety of single-step, multistep, and multistage methods. A concrete example of Crank--Nicolson integrators coupled in a multirate fashion is also included to illustrate the ideas. We conclude with some discussion of future developments, such as different interface conditions and partitioned methods. [ABSTRACT FROM AUTHOR]

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 the 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 nonsmooth integrands that are not well approximated by polynomials. [ABSTRACT FROM AUTHOR]

The numerical integration over a planar domain that is cut by an implicitly defined boundary curve is an important problem that arises, for example, in unfitted finite element methods and in isogeometric analysis on trimmed computational domains. In this paper, we introduce a a very general version of the transport theorem for moving domains defined by implicitly defined curves and use it to establish an efficient and accurate quadrature rule for this class of domains. In numerical experiments it is shown that the method achieves high orders of convergence. Our approach is suited for high-order geometrically unfitted finite element methods as well as for high-order trimmed isogeometric analysis. [ABSTRACT FROM AUTHOR]

We present a second-order accurate numerical method for a class of nonlocal nonlinear conservation laws called the "nonlocal pair-interaction model," which was recently introduced by Du, Huang, and LeFloch [SIAM J. Numer. Anal., 55 (2017), pp. 2465-2489]. Our numerical method uses second-order accurate reconstruction-based schemes for local conservation laws in conjunction with appropriate numerical integration. We show that the resulting method is total variation diminishing (TVD) and converges towards a weak solution. In fact, in contrast to local conservation laws, our second-order reconstruction-based method converges towards the unique entropy solution provided that the nonlocal interaction kernel satisfies a certain growth condition near zero. Furthermore, as the nonlocal horizon parameter in our method approaches zero we recover a well-known second-order method for local conservation laws. In addition, we answer several questions from the paper by Du, Huang, and LeFloch [SIAM J. Numer. Anal., 55 (2017), pp. 2465-2489] concerning regularity of solutions. In particular, we prove that any discontinuity present in a weak solution must be stationary and that if the interaction kernel satisfies a certain growth condition, then weak solutions are unique. We present a series of numerical experiments in which we investigate the accuracy of our second-order scheme, demonstrate shock formation in the nonlocal pair-interaction model, and examine how the regularity of the solution depends on the choice of flux function. [ABSTRACT FROM AUTHOR]

We consider the numerical approximation of Maxwell's equations in the time domain by a second-order H(curl) conforming finite element approximation. In order to enable the efficient application of explicit time-stepping schemes, we utilize a mass-lumping strategy resulting from numerical integration in conjunction with the finite element spaces introduced in [A. Elmkies and P. Joly, C. R. Acad. Sci. Paris Sér. I Math., 325 (1997), pp. 1217-1222]. We prove that this method is second-order accurate if the true solution is divergence free but the order of accuracy reduces to one in the general case. We then propose a modification of the finite element space, which yields second-order accuracy in the general case. [ABSTRACT FROM AUTHOR]

Numerical integration is encountered in all fields of numerical analysis and the engineering sciences. By now, various efficient and accurate quadrature rules are known; for instance, Gauss-type quadrature rules. In many applications, however, it might be impractical---if not even impossible---to obtain data to fit known quadrature rules. Often, experimental measurements are performed at equidistant or even scattered points in space or time. In this work, we propose stable high order quadrature rules for experimental data, which can accurately handle general weight functions. [ABSTRACT FROM AUTHOR]

We consider the numerical approximation of single-phase flow in porous media by a mixed finite element method with mass lumping. Our work extends previous results of Wheeler and Yotov, who showed that inexact numerical integration together with an appropriate choice of basis leads to some sort of mass lumping which allows eliminating the flux variables locally and reducing the mixed problem to a finite volume discretization for the pressure only. Here we construct secondorder approximations for hybrid meshes in two and three space dimensions which, similar to the method of Wheeler and Yotov, allow the local elimination of the flux variables. A full convergence analysis of the method is given for which new arguments and, in part, also new quadrature rules and finite elements are required. The analysis formally also covers approximations of higher order, and third-order elements for triangles and quadrilaterals are presented. [ABSTRACT FROM AUTHOR]

We propose a new fictitious domain finite element method, well suited for elliptic problems posed in a domain given by a level-set function without requiring a mesh fitting the boundary. To impose the Dirichlet boundary conditions, we search the approximation to the solution as a product of a finite element function with the given level-set function, which is also approximated by finite elements. Unlike other recent fictitious domain-type methods (XFEM, CutFEM), our approach does not need any nonstandard numerical integration (on cut mesh elements or on the actual boundary). We consider the Poisson equation discretized with piecewise polynomial Lagrange finite elements of any order and prove the optimal convergence of our method in the H¹ norm. Moreover, the discrete problem is proven to be well conditioned, i.e., the condition number of the associated finite element matrix is of the same order as that of a standard finite element method on a comparable conforming mesh. Numerical results confirm the optimal convergence in both H¹ and L² norms. [ABSTRACT FROM AUTHOR]

The value of a highly oscillatory integral is typically determined asymptotically by the behavior of the integrand near a small number of critical points. These include the endpoints of the integration domain and the so-called stationary points or saddle points--roots of the derivative of the phase of the integrand--where the integrand is locally nonoscillatory. Modern methods for highly oscillatory quadrature exhibit numerical issues when two such saddle points coalesce. On the other hand, integrals with coalescing saddle points are a classical topic in asymptotic analysis, where they give rise to uniform asymptotic expansions in terms of the Airy function. In this paper we construct Gaussian quadrature rules that remain uniformly accurate when two saddle points coalesce. These rules are based on orthogonal polynomials in the complex plane. We analyze these polynomials, prove their existence for even degrees, and describe an accurate and efficient numerical scheme for the evaluation of oscillatory integrals with coalescing saddle points. [ABSTRACT FROM AUTHOR]

We study multivariate numerical integration of smooth functions in weighted Sobolev spaces with dominating mixed smoothness a ≥ 2 defined over the s-dimensional unit cube. We propose a new quasi-Monte Carlo (QMC)-based quadrature rule, named the extrapolated polynomial lattice rule, which achieves the almost optimal rate of convergence. Extrapolated polynomial lattice rules are constructed in two steps: (i) construction of classical polynomial lattice rules over Fb with a consecutive sizes of nodes, bm-α+1, . . ., bm, and (ii) recursive application of Richardson extrapolation to a chain of α approximate values of the integral obtained by consecutive polynomial lattice rules. We prove the existence of good extrapolated polynomial lattice rules achieving the almost optimal order of convergence of the worst-case error in Sobolev spaces with general weights. Then, by restricting to product weights, we show that such good extrapolated polynomial lattice rules can be constructed by the fast component-by-component algorithm under a computable quality criterion. The required total construction cost is of order (s+α)N logN, which improves the currently known result for the interlaced polynomial lattice rule, which is of order saN logN. We also study the dependence of the worst-case error bound on the dimension. A big advantage of our method compared to interlaced polynomial lattice rules is that the fast QMC matrix vector method can be used in this setting, while still achieving the same rate of convergence. Such a method was previously not known. Numerical experiments for test integrands support our theoretical result. [ABSTRACT FROM AUTHOR]

In this article, we consider elliptic diusion problems with an anisotropic random diffusion coefficient. We model the notable direction in terms of a random vector field and derive regularity results for the solution's dependence on the random parameter. It turns out that the decay of the vector field's Karhunen-Loève expansion entirely determines this regularity. The obtained results allow for sophisticated quadrature methods, such as the quasi-Monte Carlo method or the anisotropic sparse grid quadrature, in order to approximate quantities of interest, like the solution's mean or the variance. Numerical examples in three spatial dimensions are provided to supplement the presented theory. [ABSTRACT FROM AUTHOR]

In this paper we consider the full discretization of linear Maxwell's equations on spatial grids which are locally refined. For such problems, explicit time integration schemes become inefficient because the smallest mesh width results in a strict CFL condition. Recently locally implicit time integration methods have become popular in overcoming the problem of so-called grid-induced stiffness. Various such schemes have been proposed in the literature and have been shown to be very efficient. However, a rigorous analysis of such methods is still lacking. In fact, the available literature focuses on error bounds which are valid on a fixed spatial mesh only but deteriorate in the limit where the smallest spatial mesh size tends to zero. Moreover, some important questions cannot be answered without such an analysis. For example, there has been no study of which elements of the spatial mesh enter the CFL condition. In this paper we provide such a rigorous analysis for a locally implicit scheme proposed by Verwer [BIT, 51 (2011), pp. 427--445] based on a variational formulation and energy techniques. [ABSTRACT FROM AUTHOR]

Strong stability preserving (SSP) methods are designed primarily for time integration of nonlinear hyperbolic PDEs, for which the permissible SSP step size varies from one step to the next. We develop the first SSP linear multistep methods (of order two and three) with variable step size, and prove their optimality, stability, and convergence. The choice of step size for multistep SSP methods is an interesting problem because the allowable step size depends on the SSP coefficient, which in turn depends on the chosen step sizes. The description of the methods includes an optimal step-size strategy. We prove sharp upper bounds on the allowable step size for explicit SSP linear multistep methods and show the existence of methods with arbitrarily high order of accuracy. The effectiveness of the methods is demonstrated through numerical examples. [ABSTRACT FROM AUTHOR]

In a recent paper by the authors, it is shown that there exists a quasi--Monte Carlo (QMC) rule which achieves the best possible rate of convergence for numerical integration in a reproducing kernel Hilbert space consisting of smooth functions. In this paper we provide an explicit construction of such an optimal order QMC rule. Our approach is to exploit both the decay and the sparsity of the Walsh coefficients of the reproducing kernel simultaneously. This can be done by applying digit interlacing composition due to Dick to digital nets with large minimum Hamming and Niederreiter--Rosenbloom--Tsfasman metrics due to Chen and Skriganov. To our best knowledge, our construction gives the first QMC rule which achieves the best possible convergence in this function space. [ABSTRACT FROM AUTHOR]

Convolving the output of discontinuous Galerkin computations with symmetric smoothness-increasing accuracy-conserving (SIAC) filters can improve both smoothness and accuracy. To extend convolution to the boundaries, several one-sided spline filters have recently been developed. This paper interprets these filters as instances of a general class of position-dependent (PSIAC) spline filters that can have nonuniform knot sequences and skip B-splines of the sequence. PSIAC filters with rational knot sequences have rational coefficients. For prototype knot sequences, such as integer sequences that may have repeated entries, PSIAC filters can be expressed in symbolic form. Based on the insight that filters for shifted or scaled knot sequences are easily derived by nonuniform scaling of one prototype filter, a single filter can be reused in different locations and at different scales. Computing a value of the convolution then simplifies to forming a scalar product of a short vector with the local output data. Restating one-sided filters in this form improves both stability and efficiency compared to their original formulation via numerical integration. PSIAC filtering is demonstrated for several established and one new boundary filter. [ABSTRACT FROM AUTHOR]

We prove that a class of monotone, W1-contractive schemes for scalar conservation laws converge at a rate of Δx2 in the Wasserstein distance (W1-distance), whenever the initial data are decreasing and consist of a finite number of piecewise constants. It is shown that the Lax- Friedrichs, Enquist-Osher, and Godunov schemes are W1-contractive. Numerical experiments are presented to illustrate the main result. To the best of our knowledge, this is the first proof of secondorder convergence of any numerical method for discontinuous solutions of nonlinear conservation laws. [ABSTRACT FROM AUTHOR]

This work embeds a multilevel Monte Carlo sampling strategy into the Monte Carlo step of the ensemble Kalman filter (EnKF) in the setting of finite dimensional signal evolution and noisy discrete-time observations. The signal dynamics is assumed to be governed by a stochastic differential equation (SDE), and a hierarchy of time grids is introduced for multilevel numerical integration of that SDE. The resulting multilevel EnKF is proved to asymptotically outperform EnKF in terms of computational cost versus approximation accuracy. The theoretical results are illustrated numerically. [ABSTRACT FROM AUTHOR]

In this paper, we study optimal estimates and comparison of the coefficients in the Gegenbauer series expansion. We propose aalternative derivation of the contour integral representation of the Gegenbauer expansion coefficients which was recently derived by Cantero and Iserles [SIAM J. Numer. Anal., 50 (2012), pp. 307--327]. With this representation, we show that optimal estimates for the Gegenbauer expansion coefficients can be derived, which in particular includes Legendre coefficients as a special case. Further, we apply these estimates to establish some rigorous and computable bounds for the truncated Gegenbauer series. In addition, we compare the decay rates of the Chebyshev and Legendre coefficients. For functions whose singularity is outside or at the endpoints of the expansion interval, asymptotic behavior of the ratio of the $n$th Legendre coefficient to the $n$th Chebyshev coefficient is given, which provides us an illuminating insight for the comparison of spectral methods based on Legendre and Chebyshev expansions. [ABSTRACT FROM AUTHOR]

We prove upper bounds on the order of convergence of Frolov's cubature formula for numerical integration in function spaces of dominating mixed smoothness on the unit cube with homogeneous boundary condition. More precisely, we study worst-case integration errors for Besov B8p,θ, and Triebel-Lizorkin spaces F8p,2 and our results treat the whole range of admissible parameters (s ⩾ 1/p). In particular, we obtain upper bounds for the difficult the case of small smoothness which is given for Triebel-Lizorkin spaces F8p,2 in case 1 < θ < p < ∝ with 1/p < s ⩽ 1/θ. The presented upper bounds on the worst-case error show a completely different behavior compared to "large" smoothness s > 1/θ. In the latter case the presented upper bounds are optimal, i.e., they cannot be improved by any other cubature formula. The optimality for "small" smoothness is open. [ABSTRACT FROM AUTHOR]

This paper studies the convergence rate of randomized quasi--Monte Carlo (RQMC) for discontinuous functions, which are often of infinite variation in the sense of Hardy and Krause. It was previously known that the root mean square error (RMSE) of RQMC is only o(n-1/2) for discontinuous functions. For certain discontinuous functions in d dimensions, we prove that the RMSE of RQMC is O(n-1/2-1/(4d-2)+∈) for any ∈ > 0 and arbitrary n. If some discontinuity boundaries are parallel to some coordinate axes, the rate can be improved to O(n...) where du denotes the so-called irregular dimension, that is, the number of axes which are not parallel to the discontinuity boundaries. Moreover, this paper shows that the RMSE is O(n-1/2-1/(2d)) for certain indicator functions. [ABSTRACT FROM AUTHOR]

In this paper, we perform the numerical analysis of a new moment/generalized polynomial chaos (gPC) based approximation method under finite numerical integration. The paper addresses the impact of this constraint on the method, in particular analyzing the interplay between aliasing and truncation errors, depending on the type of functional to be represented. We set a new theoretical result defining conditions under which the iterative procedure ensures a gain after each step, putting forward the existence of a balance between aliasing (i.e., integration accuracy) and truncation errors. We demonstrate that the iterative process is viable in this context. We emphasize the existence of two regimes, an ideal one and another one for which we suggest alternatives. The method is applied to uncertainty propagation problems, i.e., here, nonlinear mappings of a single- (or multiple-)input random variables to a single-output random variable. The originality of the approach is to automatically and iteratively adapt the stochastic approximation space in order to build an accurate recursive representation of the solution. [ABSTRACT FROM AUTHOR]

In this paper, we consider the infinite-dimensional integration problem on weighted reproducing kernel Hilbert spaces with norms induced by an underlying function space decomposition of analysis of variance type. The weights model the relative importance of different groups of variables. We present new randomized multilevel algorithms to tackle this integration problem and prove upper bounds for their randomized error. Furthermore, we provide in this setting the first non-trivial lower error bounds for general randomized algorithms, which, in particular, may be adaptive or nonlinear. These lower bounds show that our multilevel algorithms are optimal. Our analysis refines and extends the analysis provided in [F. J. Hickernell, T. Müller-Gronbach, B. Niu, and K. Ritter, J. Complexity, 26 (2010), pp. 229-254], and our error bounds improve substantially on the error bounds presented there. As an illustrative example, we discuss the unanchored Sobolev space and employ randomized quasi--Monte Carlo multilevel algorithms based on scrambled polynomial lattice rules. [ABSTRACT FROM AUTHOR]

We study modified trigonometric integrators, which generalize the popular class of trigonometric integrators for highly oscillatory Hamiltonian systems by allowing the fast frequencies to be modified. Among all methods of this class, we show that the IMEX (implicit-explicit) method, which is equivalent to applying the midpoint rule to the fast, linear part of the system and the leapfrog (Störmer/Verlet) method to the slow, nonlinear part, is distinguished by the following properties: (i) it is symplectic; (ii) it is free of artificial resonances; (iii) it is the unique method that correctly captures slow energy exchange to leading order; (iv) it conserves the total energy and a modified oscillatory energy up to to second order; (v) it is uniformly second-order accurate in the slow components; and (vi) it has the correct magnitude of deviations of the fast oscillatory energy, which is an adiabatic invariant. These theoretical results are supported by numerical experiments on the Fermi--Pasta--Ulam problem and indicate that the IMEX method, for these six properties, dominates the class of modified trigonometric integrators. [ABSTRACT FROM AUTHOR]

We present an age-structured population model that accounts for complex life cycles and competition for resources limiting the transition to maturity. Taking all of these into account leads to a new mathematical model with a state-dependent delay that cannot be analyzed using the established framework of functional differential equations or directly simulated by standard numerical schemes for age-structured populations. In this paper we present the derivation of the model and a numerical scheme to integrate the equations. Convergence of the method is proven by verifying first consistency and then stability. The difficulties involved in this proof, as well as in the implementation of the scheme, stem from the state-dependent delay term. The numerical scheme is shown to be of order 2 in the given norm, and at least of order 3/2 in the supremum norm over the mesh point values. [ABSTRACT FROM AUTHOR]

The existence of (standard) Gauss quadrature rules with respect to a nonnegative measure dµ, with support on the real axis easily can be shown with the aid of orthogonal polynomials with respect to this measure. Efficient algorithms for computing the nodes and weights of an n-point Gauss rule use the n x n symmetric tridiagonal matrix determined by the recursion coefficients for the first n orthonormal polynomials. Many rational functions that are orthogonal with respect to the measure dµ, and have real or complex conjugate poles also satisfy a short recursion relations. This paper describes how banded matrices determined by the recursion coefficients for these orthonormal rational functions can be used to efficiently compute the nodes and weights of rational Gauss quadrature rules. [ABSTRACT FROM AUTHOR]

We study the error of Gauss--Turán quadrature formulae when functions which are analytic on a neighborhood of the set of integration are considered. A computable upper bound of the error is presented which is valid for arbitrary weight functions. A comparison is made with the exact error, and a number of numerical examples for arbitrary weight functions are given which show the advantages of using such rules as well as the sharpness of the error bound. Also, a comparison is made with the other effective error bounds for some special weight functions appearing in the literature. Asymptotic error estimates when the number of nodes in the quadrature increases are presented. A couple of numerical examples which show their efficiency are included. [ABSTRACT FROM AUTHOR]

A fully discrete approximation of the linear stochastic wave equation driven by additive noise is presented. A standard finite element method is used for the spatial discretization and a stochastic trigonometric scheme for the temporal approximation. This explicit time integrator allows for error bounds independent of the space discretization and thus does not have a step-size restriction as in the often used Stormer--Verlet-leap-frog scheme. Moreover, it enjoys a trace formula as does the exact solution of our problem. These favorable properties are demonstrated with numerical experiments. [ABSTRACT FROM AUTHOR]

We consider the problem of numerical integration of a function over a discrete surface to high-order accuracy. Surface integration is a common operation in numerical computations for scientific and engineering problems. Integration over discrete surfaces (such as a surface triangula-tion) is typically limited to first- or second-order accuracy due to the piecewise linear approximations of the surface and the function. We present a novel method that can achieve third- and higher-order accuracy for integration over discrete surfaces. Our method combines a stabilized least squares approximation, a blending procedure based on linear shape functions, and high-degree quadrature rules. We present theoretical analysis of the accuracy of our method as well as experimental results of up to sixth-order accuracy with our method. [ABSTRACT FROM AUTHOR]

In this paper we consider the existence problem of minimal cubature formulas of degree 4k+1 for spherically symmetric integrals. We prove that, for a minimal formula in sufficiently high-dimensional space, there exists some concentric sphere on which the inner product of any two distinct points is rational. By using this result, we prove that for any d ≥ 2 there exist no d-dimensional minimal formulas of degrees 13 and 21 for some special integral. [ABSTRACT FROM AUTHOR]

Galerkin discretizations of integral operators in Rd require the evaluation of integrals Due to image rights restrictions, multiple line equation(s) cannot be graphically displayed., where S(1), S(1)) are d-dimensional simplices and f has a singularity at x = y. In [A. Chernov, T. von Petersdorff, and C. Schwab, M2AN Math. Model. Numer. Anal., 45 (2011), pp. 387-422] we constructed a family of hp-quadrature rules QN with N function evaluations for a class of integrands f allowing for algebraic singularities at x = y, possibly nonintegrable with respect to either dx or dy (hypersingular kernels) and Gevrey-δ smooth for x ≠ y. This is satisfied for kernels from broad classes of pseudodifferential operators. We proved that QN achieves the exponential convergence rate O(exp(-rNγ)) with the exponent γ = 1/(2dδ + 1). In this paper we consider a special singularity ||x - y||α with real α which appears frequently in appplication and prove that an improved convergence rate with γ = 1/(2dδ) is achieved if a certain one-dimensional Gauss-Jacobi quadrature rule is used in the (univariate) "singular coordinate." We also analyze approximation by tensor Gauss--Jacobi quadratures in the "regular coordinates." We illustrate the performance of the new Gauss--Jacobi rules on several numerical examples and compare it to the hp-quadratures from [A. Chernov, T. von Petersdorff, and C. Schwab, M2AN Math. Model. Numer. Anal., 45 (2011), pp. 387-422]. [ABSTRACT FROM AUTHOR]

Quasi-Monte Carlo rules are equal weight quadrature rules defined over the domain [0, 1]s. Here we introduce quasi-Monte Carlo-type rules for numerical integration of functions defined on ℝs. These rules are obtained by way of some transformation of digital nets such that locally one obtains quasi Monte Carlo rules, but at the same time, globally one also has the required distribution. We prove that these rules are optimal for numerical integration in spaces of bounded fractional variation. The analysis is based on certain tilings of the Walsh phase plane. As a proof of concept, some numerical examples are included for dimensions between 3 and 10. We compare our method with quasi-Monte Carlo rules transformed to ℝs using the inverse cumulative distribution function. In these examples, the new method significantly improves upon both methods for dimensions up to 5; no improvement is seen for dimension 10. [ABSTRACT FROM AUTHOR]

Classical Gaussian quadrature is anchored in the space of real polynomials at polynomial degree 0 and is monotonically directed by successively increasing polynomial degree. The results of recent research, investigating weighing anchor and charting various courses in search of the best currents in the space of real Laurent polynomials, are presented. Standard error bound minimization anchoring and steering algorithms which utilize the moments of positive definite strong moment functionals to guide the construction of ordered orthogonal Laurent polynomial sequences are introduced. [ABSTRACT FROM AUTHOR]

A recent result in [B. T. Helenbrook, SIAM J. Numer. Anal., 47 (2009), pp. 1304-1318] on the nonexistence of Gauss-Lobatto cubature rules on the triangle is strengthened by establishing a lower bound for the number of nodes of such rules. A method of constructing Lobatto-type cubature rules on the triangle is given and used to construct several examples. [ABSTRACT FROM AUTHOR]

A set χN of N points on the unit sphere is a spherical t-design if the average value of any polynomial of degree at most t over χN is equal to the average value of the polynomial over the sphere. This paper considers the characterization and computation of spherical t-designs on the unit sphere 핊² ⊂ ℝ³ when N ≥ (t + 1)², the dimension of the space ℙt of spherical polynomials of degree at most t. We show how to construct well conditioned spherical designs with N ≥ (t + 1)² points by maximizing the determinant of a matrix while satisfying a system of nonlinear constraints. Interval methods are then used to prove the existence of a true spherical t-design very close to the calculated points and to provide a guaranteed interval containing the determinant. The resulting spherical designs have good geometrical properties (separation and mesh norm). We discuss the usefulness of the points for both equal weight numerical integration and polynomial interpolation on the sphere and give an example. [ABSTRACT FROM AUTHOR]

If the nodes for the spectral element method are chosen to be the Gauss-Legendre-Lobatto points and a Lagrange basis is used, then the resulting mass matrix is diagonal and the method is sometimes then described as the Gauss-point mass lumped finite element scheme. We study the dispersive behavior of the scheme in detail and provide both a qualitative description of the nature of the dispersive and dissipative behavior of the scheme along with precise quantitative statements of the accuracy in terms of the mesh-size and the order of the scheme. We prove that (a) the Gauss-point mass lumped scheme (i.e., spectral element method) tends to exhibit phase lag whereas the (consistent) finite element scheme tends to exhibit phase lead; (b) the absolute accuracy of the spectral element scheme is 1/p times better than that of the finite element scheme despite the use of numerical integration; (c) when the order p, the mesh-size h, and the frequency of the wave ω satisfy 2p + 1 ≈ ωh the true wave is essentially fully resolved. As a consequence, one obtains a proof of the general rule of thumb sometimes quoted in the context of spectral element methods: π modes per wavelength are needed to resolve a wave. [ABSTRACT FROM AUTHOR]

We introduce a new class of exponential integrators for the numerical integration of large-scale systems of stiff differential equations. These so-called Rosenbrock-type methods linearize the flow in each time step and make use of the matrix exponential and related functions of the Jacobian. In contrast to standard integrators, the methods are fully explicit and do not require the numerical solution of linear systems. We analyze the convergence properties of these integrators in a semigroup framework of semilinear evolution equations in Banach spaces. In particular, we derive an abstract stability and convergence result for variable step sizes. This analysis further provides the required order conditions and thus allows us to construct pairs of embedded methods. We present a third-order method with two stages, and a fourth-order method with three stages, respectively. The application of the required matrix functions to vectors are computed by Krylov subspace approximations. We briefly discuss these implementation issues, and we give numerical examples that demonstrate the efficiency of the new integrators. [ABSTRACT FROM AUTHOR]

We construct and analyze generalized Gaussian quadrature rules for integrands with endpoint singularities or near endpoint singularities. The rules have quadrature points inside the interval of integration, and the weights are all strictly positive. Such rules date back to the study of Chebyshev sets, but their use in applications has only recently been appreciated. We provide error estimates, and we show that the convergence rate is unaffected by the singularity of the integrand. We characterize the quadrature rules in terms of two families of functions that share many properties with orthogonal polynomials but that are orthogonal with respect to a discrete scalar product that, in most cases, is not known a priori. [ABSTRACT FROM AUTHOR]

Several problems of trigonometric approximation on a hexagon and a triangle are studied using the discrete Fourier transform and orthogonal polynomials of two variables. A discrete Fourier analysis on the regular hexagon is developed in detail, from which the analysis on the triangle is deduced. The results include cubature formulas and interpolation on these domains. In particular, a trigonometric Lagrange interpolation on a triangle is shown to satisfy an explicit compact formula, which is equivalent to the polynomial interpolation on a planar region bounded by Steiner's hypocycloid. The Lebesgue constant of the interpolation is shown to be in the order of (log n)2. Furthermore, a Gauss cubature is established on the hypocycloid. [ABSTRACT FROM AUTHOR]