Integral Approximation Research Articles

Optimal numerical integration and approximation of functions on ℝd equipped with Gaussian measure

Abstract We investigate the numerical approximation of integrals over $\mathbb{R}^{d}$ equipped with the standard Gaussian measure $\gamma $ for integrands belonging to the Gaussian-weighted Sobolev spaces $W^{\alpha }_{p}(\mathbb{R}^{d}, \gamma )$ of mixed smoothness $\alpha \in \mathbb{N}$ for $1 < p < \infty $. We prove the asymptotic order of the convergence of optimal quadratures based on $n$ integration nodes and propose a novel method for constructing asymptotically optimal quadratures. As for related problems, we establish by a similar technique the asymptotic order of the linear, Kolmogorov and sampling $n$-widths in the Gaussian-weighted space $L_{q}(\mathbb{R}^{d}, \gamma )$ of the unit ball of $W^{\alpha }_{p}(\mathbb{R}^{d}, \gamma )$ for $1 \leq q < p < \infty $ and $q=p=2$.

Short Series Superconductor Density of States Integral Approximations

Accurate and fast Fermi-Dirac integral approximations are used in semiconductor device simulators to compute carrier concentrations where Boltzmann statistics cannot be applied. A similar integral can be used to compute carrier densities in the subgap region of a superconductor. The primary difference between the two integrals is the density of states (DOS) functional used. Electrons and holes in semiconductors use the spherical band approximation yielding a DOS functional proportional to the square root of state energy ( <inline-formula xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:xlink="http://www.w3.org/1999/xlink"><tex-math notation="LaTeX">$E$</tex-math></inline-formula> ). Superconductors in the subgap region can be modeled with a DOS of the form <inline-formula xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:xlink="http://www.w3.org/1999/xlink"><tex-math notation="LaTeX">$(1-\frac{E^{2}}{\Delta ^{2}})^{-0.5}$</tex-math></inline-formula> . A short series approximation using Gaussian quadrature is computed for the superconductor DOS integral. Short series approximations are also applied to the Fréchet derivatives of the integral with respect to its parameters. All short series approximations will be compared against numeric integration solutions and results in a six hundred fold reduction in integration time. A table containing the short series approximation roots and weights is given. Error plots are shown for the short series approximations at different temperatures for the niobium pair-breaking potential energy.

A node elimination algorithm for cubature of high-dimensional polytopes

Node elimination is a numerical approach for obtaining cubature rules for the approximation of multivariate integrals. Beginning with a known cubature rule, nodes are selected for elimination, and a new, more efficient rule is constructed by iteratively solving the moment equations. This paper introduces a new criterion for selecting which nodes to eliminate that is based on a linearization of the moment equation. In addition, a penalized iterative solver is introduced, that ensures that weights are positive and nodes are inside the integration domain. A strategy for constructing an initial cubature rule for various polytopes in several space dimensions is described. High efficiency rules are presented for two, three and four dimensional polytopes. The new rules are compared with rules that are obtained by combining tensor products of one dimensional quadrature rules and domain transformations, as well as with known analytically constructed cubature rules.

A Second-Order Difference Scheme for Solving a Class of Fractional Differential Equations

Introduction. Increasing accuracy in the approximation of fractional integrals, as is known, is one of the urgent tasks of computational mathematics. The purpose of this study is to create and apply a second-order difference analog to approximate the fractional Riemann-Liouville integral. Its application is investigated in solving some classes of fractional differential equations. The difference analog is designed to approximate the fractional integral with high accuracy.Materials and Methods. The paper considers a second-order difference analogue for approximating the fractional Riemann-Liouville integral, as well as a class of fractional differential equations, which contains a fractional Caputo derivative in time of the order belonging to the interval (1, 2).Results. To solve the above equations, the original fractional differential equations have been transformed into a new model that includes the Riemann-Liouville fractional integral. This transformation makes it possible to solve problems efficiently using appropriate numerical methods. Then the proposed difference analogue of the second order approximation is applied to solve the transformed model problem.Discussion and Conclusions. The stability of the proposed difference scheme is proved. An a priori estimate is obtained for the problem under consideration, which establishes the uniqueness and continuous dependence of the solution on the input data. To evaluate the accuracy of the scheme and verify the experimental order of convergence, calculations for the test problem were carried out.

Gradient-based adaptive importance samplers

Importance sampling (IS) is a powerful Monte Carlo methodology for the approximation of intractable integrals, very often involving a target probability distribution. The performance of IS heavily depends on the appropriate selection of the proposal distributions where the samples are simulated from. In this paper, we propose an adaptive importance sampler, called GRAMIS, that iteratively improves the set of proposals. The algorithm exploits geometric information of the target to adapt the location and scale parameters of those proposals. Moreover, in order to allow for a cooperative adaptation, a repulsion term is introduced that favors a coordinated exploration of the state space. This translates into a more diverse exploration and a better approximation of the target via the mixture of proposals. Moreover, we provide a theoretical justification of the repulsion term. We show the good performance of GRAMIS in two problems where the target has a challenging shape and cannot be easily approximated by a standard uni-modal proposal.

On Sharp Rate of Convergence for Discretization of Integrals Driven by Fractional Brownian Motions and Related Processes with Discontinuous Integrands

We consider equidistant approximations of stochastic integrals driven by Hölder continuous Gaussian processes of order H>12\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$H>\\frac{1}{2}$$\\end{document} with discontinuous integrands involving bounded variation functions. We give exact rate of convergence in the L1\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$L^1$$\\end{document}-distance and provide examples with different drivers. It turns out that the exact rate of convergence is proportional to n1-2H\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$n^{1-2H}$$\\end{document}, which is twice as good as the best known results in the case of discontinuous integrands and corresponds to the known rate in the case of smooth integrands. The novelty of our approach is that, instead of using multiplicative estimates for the integrals involved, we apply a change of variables formula together with some facts on convex functions allowing us to compute expectations explicitly.

Potts-Cox survival regression

A Bayesian semi-parametric survival regression model with latent partitions is introduced. Its goal is to predict survival and to cluster survival patients within the context of building prognosis systems. In order to drive cluster formation on individuals, the Potts random partition model is chosen as a prior on the covariates space. For any given partition, the proposed model assumes an interval-wise exponential distribution for the baseline hazard rate. The number of intervals is unknown. It can be estimated with a fused-lasso type penalty given by a sequential double exponential prior. Estimation and inference are done with the aid of Markov chain Monte Carlo. To simplify the computations, the Laplace's integral approximation method is used to estimate some constants and to propose parameter updates within Markov chain Monte Carlo. The methodology is illustrated with an application to cancer survival.

The Half-Space Matching method for elastodynamic scattering problems in unbounded domains

In this paper, the Half-Space Matching (HSM) method, first introduced for scalar problems, is extended to elastodynamics, to solve time-harmonic 2D scattering problems, in locally perturbed infinite anisotropic homogeneous media. The HSM formulation couples a variational formulation around the perturbations with Fourier integral representations of the outgoing solution in four overlapping half-spaces. These integral representations involve outgoing plane waves, selected according to their group velocity, and evanescent waves. Numerically, the HSM method consists in a finite element discretization of the HSM formulation, together with an approximation of the Fourier integrals. Numerical results, validating the method, are presented for different materials, isotropic and anisotropic. Comparisons with the Perfectly Matched Layers (PML) method are performed for several anisotropic materials. These results highlight the robustness of the HSM method compared to the sensitivity of the PML method with respect to its parameters.

Approximation of integration over finite groups, difference sets and association schemes

Approximation of integration over finite groups, difference sets and association schemes

Hybrid functionals for periodic systems in the density functional tight-binding method

Screened range-separated hybrid (SRSH) functionals within generalized Kohn-Sham density functional theory (GKS-DFT) have been shown to restore a general $1/(r\ensuremath{\varepsilon})$ asymptotic decay of the electrostatic interaction in dielectric environments. Major achievements of SRSH include an improved description of optical properties of solids and correct prediction of polarization-induced fundamental gap renormalization in molecular crystals. The density functional tight-binding method (DFTB) is an approximate DFT that bridges the gap between first-principles methods and empirical electronic structure schemes. While purely long-range corrected RSH are already accessible within DFTB for molecular systems, this work generalizes the theoretical foundation to also include screened range-separated hybrids, with conventional pure hybrid functionals as a special case. The presented formulation and implementation is also valid for periodic boundary conditions (PBC) beyond the $\mathrm{\ensuremath{\Gamma}}$ point. To treat periodic Fock exchange and its integrable singularity in reciprocal space, we resort to techniques successfully employed by DFT, in particular a truncated Coulomb operator and the minimum image convention. Starting from the first-principles Hartree-Fock operator, we derive suitable expressions for the DFTB method, using standard integral approximations and their efficient implementation in the $\mathrm{dftb}+$ software package. Convergence behavior is investigated and demonstrated for the polyacene series as well as two- and three-dimensional materials. Benzene and pentacene molecular and crystalline systems show the correct polarization-induced gap renormalization by SRSH-DFTB at heavily reduced computational cost compared to first-principles methods.

Computational Approach for Differential Equations with Local and Nonlocal Fractional-Order Differential Operators

Laplace transform has been used for solving differential equations of fractional order either PDEs or ODEs. However, using the Laplace transform sometimes leads to solutions in Laplace space that are not readily invertible to the real domain by analytical techniques. Therefore, numerical inversion techniques are then used to convert the obtained solution from Laplace domain into time domain. Various famous methods for numerical inversion of Laplace transform are based on quadrature approximation of Bromwich integral. The key features are the contour deformation and the choice of the quadrature rule. In this work, the Gauss–Hermite quadrature method and the contour integration method based on the trapezoidal and midpoint rule are tested and evaluated according to the criteria of applicability to actual inversion problems, applicability to different types of fractional differential equations, numerical accuracy, computational efficiency, and ease of programming and implementation. The performance and efficiency of the methods are demonstrated with the help of figures and tables. It is observed that the proposed methods converge rapidly with optimal accuracy without any time instability.

Weighted averaged Gaussian quadrature rules for modified Chebyshev measures

This paper is concerned with the approximation of integrals of a real-valued integrand over the interval [−1,1] by Gauss quadrature. The averaged and optimal averaged quadrature rules ([13,21]) provide a convenient method for approximating the error in the Gauss quadrature. However, they are applicable to all integrands that are continuous on the interval [−1,1] only if their nodes are internal, i.e. if they belong to this interval.We discuss two approaches to determine averaged quadrature rules with nodes in [−1,1]: (i) truncating the Jacobi matrix associated with the optimal averaged rule, and (ii) weighting the optimal averaged quadrature rule. We consider Chebyshev measures of the first, second, and third kinds that are modified by a linear over linear rational factor, and discuss the internality of averaged, optimal averaged, and truncated optimal averaged quadrature rules. Moreover, we show that the weighting yields internal averaged rules if a weighting parameter is properly chosen, and we provide bounds for this parameter that guarantee internality. Finally, we illustrate that the weighted averaged rules give more accurate estimates of the quadrature error than the truncated optimal averaged rules.

A Mixed-Effects Model in Which the Parameters of the Autocorrelated Error Structure Can Differ between Individuals

Research in psychology has seen a rapid increase in the usage of experience sampling methods and daily diary methods. The data that result from using these methods are typically analyzed with a mixed-effects or a multilevel model because it allows testing hypotheses about the time course of the longitudinally assessed variable or the influence of time-varying predictors in a simple way. Here, we describe an extension of this model that does not only allow to include random effects for the mean structure but also for the residual variance, for the parameter of an autoregressive process of order 1 and/or the parameter of a moving average process of order 1. After we have introduced this extension, we show how to estimate the parameters with maximum likelihood. Because the likelihood function contains complex integrals, we suggest using adaptive Gauss-Hermite quadrature and Quasi-Monte Carlo integration to approximate it. We illustrate the models using a real data example and also report the results of a small simulation study in which the two integral approximation methods are compared.

High-order asymptotic expansions of Gaussian quadrature rules with classical and generalized weight functions

Gaussian quadrature rules are a classical tool for the numerical approximation of integrals with smooth integrands and positive weight functions. We derive and explicitly list asymptotic expressions for the points and weights of Gaussian quadrature rules for three general classes of positive weight functions: analytic functions on a bounded interval with algebraic singularities at the endpoints, analytic weight functions on the halfline with exponential decay at infinity and an algebraic singularity at the finite endpoint, and analytic functions on the real line with exponential decay in both directions at infinity. The results include the Gaussian rules of classical orthogonal polynomials (Legendre, Jacobi, Laguerre and Hermite) as special cases. Explicit expressions for these cases are included in the appendix. We present experiments indicating the range of the number of points at which these expressions achieve high precision. We provide an algorithm that can compute arbitrarily many terms in these expansions for the classical cases, and many though not all terms for the generalized cases.

Corrigendum: Error bounds and asymptotic expansions for Toeplitz product functionals of unbounded spectra

We investigate error orders for integral limit approximations to traces of products of Toeplitz matrices generated by integrable functions on having some singularities at the origin. Even though a sharp error order of the above approximation is derived in Theorem 2 of Lieberman and Phillips (2004, Journal of Time Series Analysis, 25(5) 733–753), its proof contains an inaccuracy. In the present article, we reinvestigate error orders of the above trace approximation problem and rigorously validate the sharp error order derived in Lieberman and Phillips (2004, Journal of Time Series Analysis, 25(5) 733–753).

On the integral solution of elliptic Kepler’s equation

In a recent paper, Philcox, Goodman and Slepian obtain an explicit solution of the elliptic Kepler’s equation (KE) as a quotient of two contour integrals along a Jordan curve mathcal {C} = mathcal {C}(M,e) that contains the unique real solution of KE but not includes other complex zeros of KE in its interior. The aim of this paper is to study the main issues that arise in the practical implementation of this integral solution. Thus, after a study of the complex zeros of KE, several families of Jordan contours mathcal {C} = mathcal {C}(M,e) that are suitable for this integral solution are proposed. Since contours with minimal length turn out to be the more accurate for numerical purposes, several families that minimize their length are constructed. Secondly, the approximation of the contour integrals by the composite trapezoidal rule is considered. Recall that this rule is employed in the fast Fourier transform and, in spite of its lower order, displays a spectral convergence as a function of the number of nodes, which implies a very fast convergence. Finally, the results of some numerical experiments are presented to show that such a combination of appropriate contours with the composite trapezoidal rule leads to a powerful numerical method to solve KE with any desired accuracy for all values of eccentricity.

Efficient Time-Stepping for Numerical Integration Using Reinforcement Learning

Many problems in science and engineering require an efficient numerical approximation of integrals or solutions to differential equations. For systems with rapidly changing dynamics, an equidistant discretization is often inadvisable as it results in prohibitively large errors or computational effort. To this end, adaptive schemes, such as solvers based on Runge–Kutta pairs, have been developed which adapt the step size based on local error estimations at each step. While the classical schemes apply very generally and are highly efficient on regular systems, they can behave suboptimally when an inefficient step rejection mechanism is triggered by structurally complex systems such as chaotic systems. To overcome these issues, we propose a method to tailor numerical schemes to the problem class at hand. This is achieved by combining simple, classical quadrature rules or ODE solvers with data-driven time-stepping controllers. Compared with learning solution operators to ODEs directly, it generalizes better to unseen initial data as our approach employs classical numerical schemes as base methods. At the same time it can make use of identified structures of a problem class and, therefore, outperforms state-of-the-art adaptive schemes. Several examples demonstrate superior efficiency. Source code is available at https://github.com/lueckem/quadrature-ML.

Target-aware Bayesian inference via generalized thermodynamic integration

In Bayesian inference, we are usually interested in the numerical approximation of integrals that are posterior expectations or marginal likelihoods (a.k.a., Bayesian evidence). In this paper, we focus on the computation of the posterior expectation of a function f(textbf{x}). We consider a target-aware scenario where f(textbf{x}) is known in advance and can be exploited in order to improve the estimation of the posterior expectation. In this scenario, this task can be reduced to perform several independent marginal likelihood estimation tasks. The idea of using a path of tempered posterior distributions has been widely applied in the literature for the computation of marginal likelihoods. Thermodynamic integration, path sampling and annealing importance sampling are well-known examples of algorithms belonging to this family of methods. In this work, we introduce a generalized thermodynamic integration (GTI) scheme which is able to perform a target-aware Bayesian inference, i.e., GTI can approximate the posterior expectation of a given function. Several scenarios of application of GTI are discussed and different numerical simulations are provided.

An Approximation Method for Variational Inequality with Uncertain Variables

In this paper, a Stieltjes integral approximation method for uncertain variational inequality problem (UVIP) is studied. Firstly, uncertain variables are introduced on the basis of variational inequality. Since the uncertain variables are based on nonadditive measures, there is usually no density function. Secondly, the expected value model of UVIP is established after the expected value is discretized by the Stieltjes integral. Furthermore, a gap function is constructed to transform UVIP into an uncertain constraint optimization problem, and the optimal value of the constraint problem is proved to be the solution of UVIP. Finally, the convergence of solutions of the Stieltjes integral discretization approximation problem is proved.

Gauss Quadrature for Integrals and Sums

Gauss quadrature integral approximation is extended to include integrals with a measure consisting of continuous as well as discrete components. That is, we give an approximation for the integral of a function plus its sum over a discrete weighted set.