NAMar 9, 2015
Higher order QMC Galerkin discretization for parametric operator equationsJosef Dick, Frances Y. Kuo, Quoc T. Le Gia et al.
We construct quasi-Monte Carlo methods to approximate the expected values of linear functionals of Galerkin discretizations of parametric operator equations which depend on a possibly infinite sequence of parameters. Such problems arise in the numerical solution of differential and integral equations with random field inputs. We analyze the regularity of the solutions with respect to the parameters in terms of the rate of decay of the fluctuations of the input field. If $p\in (0,1]$ denotes the "summability exponent" corresponding to the fluctuations in affine-parametric families of operators, then we prove that deterministic "interlaced polynomial lattice rules" of order $α= \lfloor 1/p \rfloor+1$ in $s$ dimensions with $N$ points can be constructed using a fast component-by-component algorithm, in $\mathcal{O}(α\,s\, N\log N + α^2\,s^2 N)$ operations, to achieve a convergence rate of $\mathcal{O}(N^{-1/p})$, with the implied constant independent of $s$. This dimension-independent convergence rate is superior to the rate $\mathcal{O}(N^{-1/p+1/2})$, for $2/3\leq p\leq 1$ recently established for randomly shifted lattice rules under comparable assumptions. In our analysis we use a non-standard Banach space setting and introduce "smoothness-driven product and order dependent (SPOD)" weights for which we show fast CBC construction.
NANov 24, 2016
Multilevel higher order Quasi-Monte Carlo Bayesian EstimationJosef Dick, Robert N. Gantner, Quoc T. Le Gia et al.
We propose and analyze deterministic multilevel approximations for Bayesian inversion of operator equations with uncertain distributed parameters, subject to additive Gaussian measurement data. The algorithms use a multilevel (ML) approach based on deterministic, higher order quasi-Monte Carlo (HoQMC) quadrature for approximating the high-dimensional expectations, which arise in the Bayesian estimators, and a Petrov-Galerkin (PG) method for approximating the solution to the underlying partial differential equation (PDE). This extends the previous single-level approach from [J. Dick, R. N. Gantner, Q. T. Le Gia and Ch. Schwab, Higher order Quasi-Monte Carlo integration for Bayesian Estimation. Report 2016-13, Seminar for Applied Mathematics, ETH Zürich (in review)]. We obtain sufficient conditions which allow us to achieve arbitrarily high, algebraic convergence rates in terms of work, which are independent of the dimension of the parameter space. The convergence rates are limited only by the spatial regularity of the forward problem,the discretization order achieved by the Petrov Galerkin discretization, and by the sparsity of the uncertainty parametrization. We provide detailed numerical experiments for linear elliptic problems in two space dimensions, with $s=1024$ parameters characterizing the uncertain input, confirming the theory and showing that the ML HoQMC algorithms outperform, in terms of error vs.~computational work, both multilevel Monte Carlo (MLMC) methods and single-level (SL) HoQMC methods.
NAFeb 24, 2016
Higher order Quasi-Monte Carlo integration for Bayesian EstimationJosef Dick, Robert N. Gantner, Quoc T. Le Gia et al.
We analyze combined Quasi-Monte Carlo quadrature and Finite Element approximations in Bayesian estimation of solutions to countably-parametric operator equations with holomorphic dependence on the parameters as considered in [Cl.~Schillings and Ch.~Schwab: Sparsity in Bayesian Inversion of Parametric Operator Equations. Inverse Problems, {\bf 30}, (2014)]. Such problems arise in numerical uncertainty quantification and in Bayesian inversion of operator equations with distributed uncertain inputs, such as uncertain coefficients, uncertain domains or uncertain source terms and boundary data. We show that the parametric Bayesian posterior densities belong to a class of weighted Bochner spaces of functions of countably many variables, with a particular structure of the QMC quadrature weights: up to a (problem-dependent, and possibly large) finite dimension $S$ product weights can be used, and beyond this dimension, weighted spaces with so-called SPOD weights are used to describe the solution regularity. We establish error bounds for higher order Quasi-Monte Carlo quadrature for the Bayesian estimation based on [J.~Dick, Q.T.~LeGia and Ch.~Schwab, Higher order Quasi-Monte Carlo integration for holomorphic, parametric operator equations, Report 2014-23, SAM, ETH Zürich]. It implies, in particular, regularity of the parametric solution and of the countably-parametric Bayesian posterior density in SPOD weighted spaces. This, in turn, implies that the Quasi-Monte Carlo quadrature methods in [J. Dick, F.Y.~Kuo, Q.T.~Le Gia, D.~Nuyens, Ch.~Schwab, Higher order QMC Galerkin discretization for parametric operator equations, SINUM (2014)] are applicable to these problem classes, with dimension-independent convergence rates $\calO(N^{-1/p})$ of $N$-point HoQMC approximated Bayesian estimates, where $0<p<1$ depends only on the sparsity class of the uncertain input in the Bayesian estimation.
NAApr 1, 2013
Walsh spaces containing smooth functions and quasi-Monte Carlo rules of arbitrary high orderJosef Dick
We define a Walsh space which contains all functions whose partial mixed derivatives up to order $δ\ge 1$ exist and have finite variation. In particular, for a suitable choice of parameters, this implies that certain Sobolev spaces are contained in these Walsh spaces. For this Walsh space we then show that quasi-Monte Carlo rules based on digital $(t,α,s)$-sequences achieve the optimal rate of convergence of the worst-case error for numerical integration. This rate of convergence is also optimal for the subspace of smooth functions. Explicit constructions of digital $(t,α,s)$-sequences are given hence providing explicit quasi-Monte Carlo rules which achieve the optimal rate of convergence of the integration error for arbitrarily smooth functions.
NANov 25, 2012
Approximation of analytic functions in Korobov spacesJosef Dick, Peter Kritzer, Friedrich Pillichshammer et al.
We study multivariate $L_2$-approximation for a weighted Korobov space of analytic periodic functions for which the Fourier coefficients decay exponentially fast. The weights are defined, in particular, in terms of two sequences $\boldsymbol{a} =\{a_j\}$ and $\boldsymbol{b} =\{b_j\}$ of numbers no less than one. Let $e^{L_2-\mathrm{app},Λ}(n,s)$ be the minimal worst-case error of all algorithms that use $n$ information functionals from the class $Λ$ in the $s$-variate case. We consider two classes $Λ$: the class $Λ^{\rm all}$ consists of all linear functionals and the class $Λ^{\rm std}$ consists of only function valuations. We study (EXP) exponential convergence. This means that $$ e^{L_2-\mathrm{app},Λ}(n,s) \le C(s)\,q^{\,(n/C_1(s))^{p(s)}}\quad{for all}\quad n, s \in \mathbb{N} $$ where $q\in(0,1)$, and $C,C_1,p:\mathbb{N} \rightarrow (0,\infty)$. If we can take $p(s)=p>0$ for all $s$ then we speak of (UEXP) uniform exponential convergence. We also study EXP and UEXP with (WT) weak, (PT) polynomial and (SPT) strong polynomial tractability. These concepts are defined as follows. Let $n(\e,s)$ be the minimal $n$ for which $e^{L_2-\mathrm{app},Λ}(n,s)\le \e$. Then WT holds iff $\lim_{s+\log\,\e^{-1}\to\infty}(\log n(\e,s))/(s+\log\,\e^{-1})=0$, PT holds iff there are $c,τ_1,τ_2$ such that $n(\e,s)\le cs^{τ_1}(1+\log\,\e^{-1})^{τ_2}$ for all $s$ and $\e\in(0,1)$, and finally SPT holds iff the last estimate holds for $τ_1=0$. The infimum of $τ_2$ for which SPT holds is called the exponent $τ^*$ of SPT. We prove that the results are the same for both classes $Λ$, and obtain conditions for WT, PT, SPT with and without EXP and UEXP.
NAJan 26, 2015
Fast QMC matrix-vector multiplicationJosef Dick, Frances Y. Kuo, Quoc T. Le Gia et al.
Quasi-Monte Carlo (QMC) rules $1/N \sum_{n=0}^{N-1} f(\boldsymbol{y}_n A)$ can be used to approximate integrals of the form $\int_{[0,1]^s} f(\boldsymbol{y} A) \,\mathrm{d} \boldsymbol{y}$, where $A$ is a matrix and $\boldsymbol{y}$ is row vector. This type of integral arises for example from the simulation of a normal distribution with a general covariance matrix, from the approximation of the expectation value of solutions of PDEs with random coefficients, or from applications from statistics. In this paper we design QMC quadrature points $\boldsymbol{y}_0, ..., \boldsymbol{y}_{N-1} \in [0,1]^s$ such that for the matrix $Y = (\boldsymbol{y}_{0}^\top, ..., \boldsymbol{y}_{N-1}^\top)^\top$ whose rows are the quadrature points, one can use the fast Fourier transform to compute the matrix-vector product $Y \boldsymbol{a}^\top$, $\boldsymbol{a} \in \mathbb{R}^s$, in $\mathcal{O}(N \log N)$ operations and at most $s-1$ extra additions. The proposed method can be applied to lattice rules, polynomial lattice rules and a certain type of Korobov $p$-set. The approach is illustrated computationally by three numerical experiments. The first test considers the generation of points with normal distribution and general covariance matrix, the second test applies QMC to high-dimensional, affine-parametric, elliptic partial differential equations with uniformly distributed random coefficients, and the third test addresses Finite-Element discretizations of elliptic partial differential equations with high-dimensional, log-normal random input data. All numerical tests show a significant speed-up of the computation times of the fast QMC matrix method compared to a conventional implementation as the dimension becomes large.
NAApr 1, 2013
Explicit constructions of quasi-Monte Carlo rules for the numerical integration of high dimensional periodic functionsJosef Dick
In this paper we give explicit constructions of point sets in the $s$ dimensional unit cube yielding quasi-Monte Carlo algorithms which achieve the optimal rate of convergence of the worst-case error for numerically integrating high dimensional periodic functions. In the classical measure $P_α$ of the worst-case error introduced by Korobov the convergence is of $\landau(N^{-\min(α,d)} (\log N)^{sα-2})$ for every even integer $α\ge 1$, where $d$ is a parameter of the construction which can be chosen arbitrarily large and $N$ is the number of quadrature points. This convergence rate is known to be best possible up to some $\log N$ factors. We prove the result for the deterministic and also a randomized setting. The construction is based on a suitable extension of digital $(t,m,s)$-nets over the finite field $\integer_b$.
NANov 16, 2012
Lattice rules for nonperiodic smooth integrandsJosef Dick, Dirk Nuyens, Friedrich Pillichshammer
The aim of this paper is to show that one can achieve convergence rates of $N^{-α+ δ}$ for $α> 1/2$ (and for $δ> 0$ arbitrarily small) for nonperiodic $α$-smooth cosine series using lattice rules without random shifting. The smoothness of the functions can be measured by the decay rate of the cosine coefficients. For a specific choice of the parameters the cosine series space coincides with the unanchored Sobolev space of smoothness 1. We study the embeddings of various reproducing kernel Hilbert spaces and numerical integration in the cosine series function space and show that by applying the so-called tent transformation to a lattice rule one can achieve the (almost) optimal rate of convergence of the integration error. The same holds true for symmetrized lattice rules for the tensor product of the direct sum of the Korobov space and cosine series space, but with a stronger dependence on the dimension in this case.
NAJan 24, 2011
A simple Proof of Stolarsky's Invariance PrincipleJohann S. Brauchart, Josef Dick
Stolarsky [Proc. Amer. Math. Soc. 41 (1973), 575--582] showed a beautiful relation that balances the sums of distances of points on the unit sphere and their spherical cap $\mathbb{L}_2$-discrepancy to give the distance integral of the uniform measure on the sphere a potential-theoretical quantity (Bj{ö}rck [Ark. Mat. 3 (1956), 255--269]). Read differently it expresses the worst-case numerical integration error for functions from the unit ball in a certain Hilbert space setting in terms of the $\mathbb{L}_2$-discrepancy and vice versa (first author and Womersley [Preprint]). In this note we give a simple proof of the invariance principle using reproducing kernel Hilbert spaces.
NANov 20, 2012
Higher order scrambled digital nets achieve the optimal rate of the root mean square error for smooth integrandsJosef Dick
We study a random sampling technique to approximate integrals $\int_{[0,1]^s}f(\mathbf{x})\,\mathrm{d}\mathbf{x}$ by averaging the function at some sampling points. We focus on cases where the integrand is smooth, which is a problem which occurs in statistics. The convergence rate of the approximation error depends on the smoothness of the function $f$ and the sampling technique. For instance, Monte Carlo (MC) sampling yields a convergence of the root mean square error (RMSE) of order $N^{-1/2}$ (where $N$ is the number of samples) for functions $f$ with finite variance. Randomized QMC (RQMC), a combination of MC and quasi-Monte Carlo (QMC), achieves a RMSE of order $N^{-3/2+\varepsilon}$ under the stronger assumption that the integrand has bounded variation. A combination of RQMC with local antithetic sampling achieves a convergence of the RMSE of order $N^{-3/2-1/s+\varepsilon}$ (where $s\ge1$ is the dimension) for functions with mixed partial derivatives up to order two. Additional smoothness of the integrand does not improve the rate of convergence of these algorithms in general. On the other hand, it is known that without additional smoothness of the integrand it is not possible to improve the convergence rate. This paper introduces a new RQMC algorithm, for which we prove that it achieves a convergence of the root mean square error (RMSE) of order $N^{-α-1/2+\varepsilon}$ provided the integrand satisfies the strong assumption that it has square integrable partial mixed derivatives up to order $α>1$ in each variable. Known lower bounds on the RMSE show that this rate of convergence cannot be improved in general for integrands with this smoothness. We provide numerical examples for which the RMSE converges approximately with order $N^{-5/2}$ and $N^{-7/2}$, in accordance with the theoretical upper bound.
NAAug 5, 2014
Construction of interlaced scrambled polynomial lattice rules of arbitrary high orderTakashi Goda, Josef Dick
Higher order scrambled digital nets are randomized quasi-Monte Carlo rules which have recently been introduced in [J. Dick, Ann. Statist., 39 (2011), 1372--1398] and shown to achieve the optimal rate of convergence of the root mean square error for numerical integration of smooth functions defined on the $s$-dimensional unit cube. The key ingredient there is a digit interlacing function applied to the components of a randomly scrambled digital net whose number of components is $ds$, where the integer $d$ is the so-called interlacing factor. In this paper, we replace the randomly scrambled digital nets by randomly scrambled polynomial lattice point sets, which allows us to obtain a better dependence on the dimension while still achieving the optimal rate of convergence. Our results apply to Owen's full scrambling scheme as well as the simplifications studied by Hickernell, Matoušek and Owen. We consider weighted function spaces with general weights, whose elements have square integrable partial mixed derivatives of order up to $α\ge 1$, and derive an upper bound on the variance of the estimator for higher order scrambled polynomial lattice rules. Employing our obtained bound as a quality criterion, we prove that the component-by-component construction can be used to obtain explicit constructions of good polynomial lattice point sets. By first constructing classical polynomial lattice point sets in base $b$ and dimension $ds$, to which we then apply the interlacing scheme of order $d$, we obtain a construction cost of the algorithm of order $\mathcal{O}(dsmb^m)$ operations using $\mathcal{O}(b^m)$ memory in case of product weights, where $b^m$ is the number of points in the polynomial lattice point set.
NAMay 18, 2011
Efficient calculation of the worst-case error and (fast) component-by-component construction of higher order polynomial lattice rulesJan Baldeaux, Josef Dick, Gunther Leobacher et al.
We show how to obtain a fast component-by-component construction algorithm for higher order polynomial lattice rules. Such rules are useful for multivariate quadrature of high-dimensional smooth functions over the unit cube as they achieve the near optimal order of convergence. The main problem addressed in this paper is to find an efficient way of computing the worst-case error. A general algorithm is presented and explicit expressions for base~2 are given. To obtain an efficient component-by-component construction algorithm we exploit the structure of the underlying cyclic group. We compare our new higher order multivariate quadrature rules to existing quadrature rules based on higher order digital nets by computing their worst-case error. These numerical results show that the higher order polynomial lattice rules improve upon the known constructions of quasi-Monte Carlo rules based on higher order digital nets.
NANov 15, 2017
On the optimal order of integration in Hermite spaces with finite smoothnessJosef Dick, Christian Irrgeher, Gunther Leobacher et al.
We study the numerical approximation of integrals over $\mathbb{R}^s$ with respect to the standard Gaussian measure for integrands which lie in certain Hermite spaces of functions. The decay rate of the associated sequence is specified by a single integer parameter which determines the smoothness classes and the inner product can be expressed via $L_2$ norms of the derivatives of the function. We map higher order digital nets from the unit cube to a suitable subcube of $\mathbb{R}^s$ via a linear transformation and show that such rules achieve, apart from powers of $\log N$, the optimal rate of convergence of the integration error.
NASep 2, 2014
Proof Techniques in Quasi-Monte Carlo TheoryJosef Dick, Aicke Hinrichs, Friedrich Pillichshammer
In this survey paper we discuss some tools and methods which are of use in quasi-Monte Carlo (QMC) theory. We group them in chapters on Numerical Analysis, Harmonic Analysis, Algebra and Number Theory, and Probability Theory. We do not provide a comprehensive survey of all tools, but focus on a few of them, including reproducing and covariance kernels, Littlewood-Paley theory, Riesz products, Minkowski's fundamental theorem, exponential sums, diophantine approximation, Hoeffding's inequality and empirical processes, as well as other tools. We illustrate the use of these methods in QMC using examples.
NAOct 16, 2012
Infinite-Dimensional Integration in Weighted Hilbert Spaces: Anchored Decompositions, Optimal Deterministic Algorithms, and Higher Order ConvergenceJosef Dick, Michael Gnewuch
We study numerical integration of functions depending on an infinite number of variables. We provide lower error bounds for general deterministic linear algorithms and provide matching upper error bounds with the help of suitable multilevel algorithms and changing dimension algorithms. More precisely, the spaces of integrands we consider are weighted reproducing kernel Hilbert spaces with norms induced by an underlying anchored function space decomposition. Here the weights model the relative importance of different groups of variables. The error criterion used is the deterministic worst case error. We study two cost models for function evaluation which depend on the number of active variables of the chosen sample points, and two classes of weights, namely product and order-dependent (POD) weights and the newly introduced weights with finite active dimension. We show for these classes of weights that multilevel algorithms achieve the optimal rate of convergence in the first cost model while changing dimension algorithms achieve the optimal convergence rate in the second model. As an illustrative example, we discuss the anchored Sobolev space with smoothness parameter $α$ and provide new optimal quasi-Monte Carlo multilevel algorithms and quasi-Monte Carlo changing dimension algorithms based on higher-order polynomial lattice rules.
NAJul 29, 2011
Quasi-Monte Carlo rules for numerical integration over the unit sphere $\mathbb{S}^2$Johann S. Brauchart, Josef Dick
We study numerical integration on the unit sphere $\mathbb{S}^2 \subset \mathbb{R}^3$ using equal weight quadrature rules, where the weights are such that constant functions are integrated exactly. The quadrature points are constructed by lifting a $(0,m,2)$-net given in the unit square $[0,1]^2$ to the sphere $\mathbb{S}^2$ by means of an area preserving map. A similar approach has previously been suggested by Cui and Freeden [SIAM J. Sci. Comput. 18 (1997), no. 2]. We prove three results. The first one is that the construction is (almost) optimal with respect to discrepancies based on spherical rectangles. Further we prove that the point set is asymptotically uniformly distributed on $\mathbb{S}^2$. And finally, we prove an upper bound on the spherical cap $L_2$-discrepancy of order $N^{-1/2} (\log N)^{1/2}$ (where $N$ denotes the number of points). This slightly improves upon the bound on the spherical cap $L_2$-discrepancy of the construction by Lubotzky, Phillips and Sarnak [Comm. Pure Appl. Math. 39 (1986), 149--186]. Numerical results suggest that the $(0,m,2)$-nets lifted to the sphere $\mathbb{S}^2$ have spherical cap $L_2$-discrepancy converging with the optimal order of $N^{-3/4}$.
NTJun 3, 2013
Optimal $\mathcal{L}_2$ discrepancy bounds for higher order digital sequences over the finite field $\mathbb{F}_2$Josef Dick, Friedrich Pillichshammer
We show that the $\mathcal{L}_2$ discrepancy of the explicitly constructed infinite sequences of points $(\boldsymbol{x}_0,\boldsymbol{x}_1, \boldsymbol{x}_2,...)$ in $[0,1)^s$ over $\mathbb{F}_2$ introduced in [J. Dick, Walsh spaces containing smooth functions and quasi-Monte Carlo rules of arbitrary high order. SIAM J. Numer. Anal., {\bf 46}, 1519--1553, 2008] satisfy $$\mathcal{L}_{2,N}(\{\boldsymbol{x}_0,\boldsymbol{x}_1,..., \boldsymbol{x}_{N-1}\}) \le C_s N^{-1} (\log N)^{s/2} \quad {for all} N \ge 2,$$ and $$\mathcal{L}_{2,2^m}(\{\boldsymbol{x}_0,\boldsymbol{x}_1,..., \boldsymbol{x}_{2^m-1}\}) \le C_s 2^{-m} m^{(s-1)/2} \quad {for all} m \ge 1,$$ where $C_s > 0$ is a constant independent of $N$ and $m$. These results are best possible by lower bounds in [P.D. Proinov, On the $L^2$ discrepancy of some infinite sequences. Serdica, {\bf 11}, 3--12, 1985] and [K. F. Roth, On irregularities of distribution. Mathematika, {\bf 1}, 73--79, 1954]. Further, for every $N \ge 2$ we explicitly construct finite point sets $\{\boldsymbol{y}_0,..., \boldsymbol{y}_{N-1}\}$ in $[0,1)^s$ such that $$\mathcal{L}_{2,N}(\{\boldsymbol{y}_0,\boldsymbol{y}_1,..., \boldsymbol{y}_{N-1}\}) \le C_s N^{-1} (\log N)^{(s-1)/2}.$$ Another solution for finite point sets by a different construction was previously shown in [W. W. L. Chen and M. M. Skriganov, Explicit constructions in the classical mean squares problem in irregularity of point distribution. J. Reine Angew. Math., {\bf 545}, 67--95, 2002].
COJan 15, 2016
Discrepancy bounds for uniformly ergodic Markov chain quasi-Monte CarloJosef Dick, Daniel Rudolf, Houying Zhu
Markov chains can be used to generate samples whose distribution approximates a given target distribution. The quality of the samples of such Markov chains can be measured by the discrepancy between the empirical distribution of the samples and the target distribution. We prove upper bounds on this discrepancy under the assumption that the Markov chain is uniformly ergodic and the driver sequence is deterministic rather than independent $U(0,1)$ random variables. In particular, we show the existence of driver sequences for which the discrepancy of the Markov chain from the target distribution with respect to certain test sets converges with (almost) the usual Monte Carlo rate of $n^{-1/2}$.
NASep 29, 2014
Higher Order Quasi Monte-Carlo Integration in Uncertainty QuantificationJosef Dick, Quoc Thong Le Gia, Christoph Schwab
We review recent results on dimension-robust higher order convergence rates of Quasi-Monte Carlo Petrov-Galerkin approximations for response functionals of infinite-dimensional, parametric operator equations which arise in computational uncertainty quantification.
NAAug 20, 2014
Spatial low-discrepancy sequences, spherical cone discrepancy, and applications in financial modelingJohann S. Brauchart, Josef Dick, Lou Fang
In this paper we introduce a reproducing kernel Hilbert space defined on $\mathbb{R}^{d+1}$ as the tensor product of a reproducing kernel defined on the unit sphere $\mathbb{S}^{d}$ in $\mathbb{R}^{d+1}$ and a reproducing kernel defined on $[0,\infty)$. We extend Stolarsky's invariance principle to this case and prove upper and lower bounds for numerical integration in the corresponding reproducing kernel Hilbert space. The idea of separating the direction from the distance from the origin can also be applied to the construction of quadrature methods. An extension of the area-preserving Lambert transform is used to generate points on $\mathbb{S}^{d-1}$ via lifting Sobol' points in $[0,1)^{d}$ to the sphere. The $d$-th component of each Sobol' point, suitably transformed, provides the distance information so that the resulting point set is normally distributed in $\mathbb{R}^{d}$. Numerical tests provide evidence of the usefulness of constructing Quasi-Monte Carlo type methods for integration in such spaces. We also test this method on examples from financial applications (option pricing problems) and compare the results with traditional methods for numerical integration in $\mathbb{R}^{d}$.
NASep 15, 2011
Point sets on the sphere $\mathbb{S}^2$ with small spherical cap discrepancyChristoph Aistleitner, Johann Brauchart, Josef Dick
In this paper we study the geometric discrepancy of explicit constructions of uniformly distributed points on the two-dimensional unit sphere. We show that the spherical cap discrepancy of random point sets, of spherical digital nets and of spherical Fibonacci lattices converges with order $N^{-1/2}$. Such point sets are therefore useful for numerical integration and other computational simulations. The proof uses an area-preserving Lambert map. A detailed analysis of the level curves and sets of the pre-images of spherical caps under this map is given.
NAOct 2, 2012
On the fast computation of the weight enumerator polynomial and the $t$ value of digital nets over finite abelian groupsJosef Dick, Makoto Matsumoto
In this paper we introduce digital nets over finite abelian groups which contain digital nets over finite fields and certain rings as a special case. We prove a MacWilliams type identity for such digital nets. This identity can be used to compute the strict $t$-value of a digital net over finite abelian groups. If the digital net has $N$ points in the $s$ dimensional unit cube $[0,1]^s$, then the $t$-value can be computed in $\mathcal{O}(N s \log N)$ operations and the weight enumerator polynomial can be computed in $\mathcal{O}(N s (\log N)^2)$ operations, where operations mean arithmetic of integers. By precomputing some values the number of operations of computing the weight enumerator polynomial can be reduced further.
NAApr 6, 2010
A Construction of Polynomial Lattice Rules with Small Gain CoefficientsJan Baldeaux, Josef Dick
In this paper we construct polynomial lattice rules which have, in some sense, small gain coefficients using a component-by-component approach. The gain coefficients, as introduced by Owen, indicate to what degree the method improves upon Monte Carlo. We show that the variance of an estimator based on a scrambled polynomial lattice rule constructed component-by-component decays at a rate of $N^{-(2α+ 1) +δ}$, for all $δ>0$, assuming that the function under consideration has bounded fractional variation of order $α$ and where $N$ denotes the number of quadrature points. An analogous result is obtained for Korobov polynomial lattice rules. It is also established that these rules are almost optimal for the function space considered in this paper. Furthermore, we discuss the implementation of the component-by-component approach and show how to reduce the computational cost associated with it. Finally, we present numerical results comparing scrambled polynomial lattice rules and scrambled digital nets.
NTApr 15, 2016
Optimal $L_p$-discrepancy bounds for second order digital sequencesJosef Dick, Aicke Hinrichs, Lev Markhasin et al.
The $L_p$-discrepancy is a quantitative measure for the irregularity of distribution modulo one of infinite sequences. In 1986 Proinov proved for all $p>1$ a lower bound for the $L_p$-discrepancy of general infinite sequences in the $d$-dimensional unit cube, but it remained an open question whether this lower bound is best possible in the order of magnitude until recently. In 2014 Dick and Pillichshammer gave a first construction of an infinite sequence whose order of $L_2$-discrepancy matches the lower bound of Proinov. Here we give a complete solution to this problem for all finite $p > 1$. We consider so-called order $2$ digital $(t,d)$-sequences over the finite field with two elements and show that such sequences achieve the optimal order of $L_p$-discrepancy simultaneously for all $p \in (1,\infty)$.
NAOct 26, 2018
Richardson extrapolation of polynomial lattice rulesJosef Dick, Takashi Goda, Takehito Yoshiki
We study multivariate numerical integration of smooth functions in weighted Sobolev spaces with dominating mixed smoothness $α\geq 2$ defined over the $s$-dimensional unit cube. We propose a new quasi-Monte Carlo (QMC)-based quadrature rule, named \emph{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 $\mathbb{F}_b$ with $α$ consecutive sizes of nodes, $b^{m-α+1},\ldots,b^{m}$, 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\log N$, which improves the currently known result for interlaced polynomial lattice rule, that is of order $sαN\log N$. 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.
NADec 21, 2015
Digital inversive vectors can achieve strong polynomial tractability for the weighted star discrepancy and for multivariate integrationJosef Dick, Domingo Gomez-Perez, Friedrich Pillichshammer et al.
We study high-dimensional numerical integration in the worst-case setting. The subject of tractability is concerned with the dependence of the worst-case integration error on the dimension. Roughly speaking, an integration problem is tractable if the worst-case error does not explode exponentially with the dimension. Many classical problems are known to be intractable. However, sometimes tractability can be shown. Often such proofs are based on randomly selected integration nodes. Of course, in applications true random numbers are not available and hence one mimics them with pseudorandom number generators. This motivates us to propose the use of pseudorandom vectors as underlying integration nodes in order to achieve tractability. In particular, we consider digital inverse vectors and present two examples of problems, the weighted star discrepancy and integration of Hölder continuous, absolute convergent Fourier- and cosine series, where the proposed method is successful.
COAug 8, 2014
Discrepancy Estimates for Acceptance-Rejection Samplers Using Stratified InputsHouying Zhu, Josef Dick
In this paper we propose an acceptance-rejection sampler using stratified inputs as diver sequence. We estimate the discrepancy of the points generated by this algorithm. First we show an upper bound on the star discrepancy of order $N^{-1/2-1/(2s)}$. Further we prove an upper bound on the $q$-th moment of the $L_q$-discrepancy $(\mathbb{E}[N^{q}L^{q}_{q,N}])^{1/q}$ for $2\le q\le \infty$, which is of order $N^{(1-1/s)(1-1/q)}$. We also present an improved convergence rate for a deterministic acceptance-rejection algorithm using $(t,m,s)-$nets as driver sequence.
NADec 14, 2016
Construction of interlaced polynomial lattice rules for infinitely differentiable functionsJosef Dick, Takashi Goda, Kosuke Suzuki et al.
We study multivariate integration over the $s$-dimensional unit cube in a weighted space of infinitely differentiable functions. It is known from a recent result by Suzuki that there exists a good quasi-Monte Carlo (QMC) rule which achieves a super-polynomial convergence of the worst-case error in this function space, and moreover, that this convergence behavior is independent of the dimension under a certain condition on the weights. In this paper we provide a constructive approach to finding a good QMC rule achieving such a dimension-independent super-polynomial convergence of the worst-case error. Specifically, we prove that interlaced polynomial lattice rules, with an interlacing factor chosen properly depending on the number of points $N$ and the weights, can be constructed using a fast component-by-component algorithm in at most $O(sN(\log N)^2)$ arithmetic operations to achieve a dimension-independent super-polynomial convergence. The key idea for the proof of the worst-case error bound is to use a variant of Jensen's inequality with a purposely-designed concave function.
NASep 23, 2011
Random weights, robust lattice rules and the geometry of the cbc$r$c algorithmJosef Dick
In this paper we study lattice rules which are cubature formulae to approximate integrands over the unit cube $[0,1]^s$ from a weighted reproducing kernel Hilbert space. We assume that the weights are independent random variables with a given mean and variance for two reasons stemming from practical applications: (i) It is usually not known in practice how to choose the weights. Thus by assuming that the weights are random variables, we obtain robust constructions (with respect to the weights) of lattice rules. This, to some extend, removes the necessity to carefully choose the weights. (ii) In practice it is convenient to use the same lattice rule for many different integrands. The best choice of weights for each integrand may vary to some degree, hence considering the weights random variables does justice to how lattice rules are used in applications. We also study a generalized version which uses $r$ constraints which we call the cbc$r$c (component-by-component with $r$ constraints) algorithm. We show that lattice rules generated by the cbc$r$c algorithm simultaneously work well for all weights in a subspace spanned by the chosen weights $\boldsymbolγ^{(1)},..., \boldsymbolγ^{(r)}$. Thus, in applications, instead of finding one set of weights, it is enough to find an $r$ dimensional convex polytope in which the optimal weights lie. The price for this method is a factor $r$ in the upper bound on the error and in the construction cost of the lattice rule. Thus the burden of determining one set of weights very precisely can be shifted to the construction of good lattice rules.
90.2NAApr 23
Spherical Cap $L_2$ Discrepancy -- Blessing of Dimensionality and a Balanced Large-Cap VariantJohann S. Brauchart, Josef Dick, Friedrich Pillichshammer
We prove that the information complexity (i.e., the inverse) of the classical spherical cap $L_2$ discrepancy on the $d$-dimensional sphere $\mathbb{S}^d$ decreases with dimension $d$, indicating a ``blessing of dimensionality'' for the associated numerical integration problem. We then introduce a modified spherical cap $L_2$ discrepancy that emphasizes large caps (close to hemispheres). For this variant, the problem does not become easier with increasing $d$. We also establish a Stolarsky invariance principle which connects the modified spherical cap $L_2$ discrepancy to numerical integration in the Sobolev space $H^{(d+1)/2}(\mathbb{S}^d)$, represented by the reproducing kernel $K(\boldsymbol{x}, \boldsymbol{y}) = 1 - \tfrac{1}{\sqrt{2}} \|\boldsymbol{x} - \boldsymbol{y}\|$. Stolarsky's invariance principle then implies that the worst-case integration error in this space grows polynomially with $d$.
NAAug 4, 2008
A Multivariate Fast Discrete Walsh Transform with an Application to Function InterpolationKwong-Ip Liu, Josef Dick, Fred J. Hickernell
For high dimensional problems, such as approximation and integration, one cannot afford to sample on a grid because of the curse of dimensionality. An attractive alternative is to sample on a low discrepancy set, such as an integration lattice or a digital net. This article introduces a multivariate fast discrete Walsh transform for data sampled on a digital net that requires only $O(N \log N)$ operations, where $N$ is the number of data points. This algorithm and its inverse are digital analogs of multivariate fast Fourier transforms. This fast discrete Walsh transform and its inverse may be used to approximate the Walsh coefficients of a function and then construct a spline interpolant of the function. This interpolant may then be used to estimate the function's effective dimension, an important concept in the theory of numerical multivariate integration. Numerical results for various functions are presented.
NAJul 16, 2014
The inverse of the star-discrepancy problem and the generation of pseudo-random numbersJosef Dick, Friedrich Pillichshammer
The inverse of the star-discrepancy problem asks for point sets $P_{N,s}$ of size $N$ in the $s$-dimensional unit cube $[0,1]^s$ whose star-discrepancy $D^\ast(P_{N,s})$ satisfies $$D^\ast(P_{N,s}) \le C \sqrt{s/N},$$ where $C> 0$ is a constant independent of $N$ and $s$. The first existence results in this direction were shown by Heinrich, Novak, Wasilkowski, and Woźniakowski in 2001, and a number of improvements have been shown since then. Until now only proofs that such point sets exist are known. Since such point sets would be useful in applications, the big open problem is to find explicit constructions of suitable point sets $P_{N,s}$. We review the current state of the art on this problem and point out some connections to pseudo-random number generators.
NAApr 15, 2016
A Discrepancy Bound for Deterministic Acceptance-Rejection Samplers Beyond $N^{-1/2}$ in Dimension 1Houying Zhu, Josef Dick
In this paper we consider an acceptance-rejection (AR) sampler based on deterministic driver sequences. We prove that the discrepancy of an $N$ element sample set generated in this way is bounded by $\mathcal{O} (N^{-2/3}\log N)$, provided that the target density is twice continuously differentiable with non-vanishing curvature and the AR sampler uses the driver sequence $$\mathcal{K}_M= \{( j α, j β) ~~ mod~~1 \mid j = 1,\ldots,M\}, $$ where $α,β$ are real algebraic numbers such that $1,α,β$ is a basis of a number field over $\mathbb{Q}$ of degree $3$. For the driver sequence $$\mathcal{F}_k= \{ ({j}/{F_k}, \{jF_{k-1}/{F_k}\} ) \mid j=1,\ldots, F_k\},$$ where $F_k$ is the $k$-th Fibonacci number and $\{x\}=x-\lfloor x \rfloor$ is the fractional part of a non-negative real number $x$, we can remove the $\log$ factor to improve the convergence rate to $\mathcal{O}(N^{-2/3})$, where again $N$ is the number of samples we accepted. We also introduce a criterion for measuring the goodness of driver sequences. The proposed approach is numerically tested by calculating the star-discrepancy of samples generated for some target densities using $\mathcal{K}_M$ and $\mathcal{F}_k$ as driver sequences. These results confirm that achieving a convergence rate beyond $N^{-1/2}$ is possible in practice using $\mathcal{K}_M$ and $\mathcal{F}_k$ as driver sequences in the acceptance-rejection sampler.
NAFeb 14, 2018
Digital net properties of a polynomial analogue of Frolov's constructionJosef Dick, Friedrich Pillichshammer, Kosuke Suzuki et al.
Frolov's cubature formula on the unit hypercube has been considered important since it attains an optimal rate of convergence for various function spaces. Its integration nodes are given by shrinking a suitable full rank $\mathbb{Z}$-lattice in $\mathbb{R}^d$ and taking all points inside the unit cube. The main drawback of these nodes is that they are hard to find computationally, especially in high dimensions.In such situations, quasi-Monte Carlo (QMC) rules based on digital nets have proven to be successful. However, there is still no construction known that leads to QMC rules which are optimal in the same generality as Frolov's. In this paper we investigate a polynomial analog of Frolov's cubature formula, which we expect to be important in this respect. This analog is defined in a field of Laurent series with coefficients in a finite field. A similar approach was previously studied in [M.~B.~Levin. Adelic constructions of low discrepancy sequences. Online Journal of Analytic Combinatorics. Issue 5, 2010.]. We show that our construction is a $(t,m,d)$-net, which also implies bounds on its star-discrepancy and the error of the corresponding cubature rule. Moreover, we show that our cubature rule is a QMC rule, whereas Frolov's is not, and provide an algorithm to determine its integration nodes explicitly. To this end we need to extend the notion of $(t,m,d)$-nets to fit the situation that the points can have infinite digit expansion and develop a duality theory. Additionally, we adapt the notion of admissible lattices to our setting and prove its significance.
NAJul 23, 2012
A fast Fourier transform method for computing the weight enumerator polynomial and trigonometric degree of lattice rulesJosef Dick
A fast Fourier transform method for computing the weight enumerator polynomial and trigonometric degree of lattice rules is introduced.
NAMar 20, 2012
A higher order Blokh-Zyablov propagation rule for higher order netsJosef Dick, Peter Kritzer
Higher order nets were introduced by Dick as a generalisation of classical $(t,m,s)$-nets, which are point sets frequently used in quasi-Monte Carlo integration algorithms. Essential tools in finding such point sets of high quality are propagation rules, which make it possible to generate new higher order nets from existing higher order nets and even classical $(t,m,s)$-nets. Such propagation rules for higher order nets were first considered by the authors in [J. Dick, P. Kritzer. Duality theory and propagation rules for generalized digital nets. Math. Comp. 79, 993--1017, 2010] and further developed in [J. Baldeaux, J. Dick, F. Pillichshammer. Duality theory and propagation rules for higher order nets. Discrete Math. 311, 362--386, 2011]. In [E.L. Blokh, V.V. Zyablov. Coding of generalized concatenated codes. Problems of Information Transmission, 10, 218--222, 1974] Blokh and Zyablov established a very general propagation rule for linear codes. This propagation rule has been extended to $(t,m,s)$-nets by Schürer and Schmid in [R. Schürer, W.Ch. Schmid. \textit{MinT---the database of optimal net, code, OA, and OOA parameters}. Available at: \texttt{http://mint.sbg.ac.at}]. In this paper we show that this propagation rule can also be extended to higher order nets. Examples indicate that this propagation rule yields new higher order nets with significantly higher quality.
NANov 11, 2010
Quasi-Monte Carlo numerical integration on $\mathbb{R}^s$: digital nets and worst-case errorJosef Dick
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 $\mathbb{R}^s$. These rules are obtained by way of some transformation of digital nets such that locally one obtains qMC rules, but at the same time, globally one also has the required distribution. We prove that these rules are optimal for numerical integration in fractional Besov type spaces. The analysis is based on certain tilings of the Walsh phase plane.
58.4NAApr 1
A deterministic multiple-shift lattice algorithm for function approximation in Korobov and half-period Cosine spacesJiarui Du, Josef Dick
Approximating multivariate periodic functions in weighted Korobov spaces via rank-1 lattices is bottlenecked by frequency aliasing. Existing optimal-rate methods rely on randomized constructions or massive pre-computations. We propose a fully deterministic multiple-shift lattice algorithm without pre pre-computations. First, we develop a uniform shift framework for aliased frequency fibers that decouples and reduces sampling costs. Second, leveraging the Chinese Remainder Theorem and the Weil bound, we introduce an adaptive hybrid construction that algebraically guarantees the full rank and bounded condition number of the reconstruction matrix. We rigorously prove that this deterministic method maintains the optimal convergence rate in the worst-case setting. Furthermore, we extend this framework to non-periodic domains via the tent transformation. By establishing a strict projection equivalence, we prove the algorithm attains optimal $L_2$ and $L_\infty$ approximation orders in the half-period cosine space, successfully resolving an open theoretical problem posed by Suryanarayana et al. (2016). This mathematically validates the proposed algorithm as a generic meshless spectral solver for high-dimensional boundary value problems, such as the Poisson equation with Neumann conditions. Numerical experiments corroborate the theoretical bounds, demonstrating an order-of-magnitude reduction in sampling complexity over probabilistic baselines while ensuring absolute deterministic stability.
NAMay 28, 2025
A decomposition-based robust training of physics-informed neural networks for nearly incompressible linear elasticityJosef Dick, Seungchan Ko, Quoc Thong Le Gia et al.
Due to divergence instability, the accuracy of low-order conforming finite element methods for nearly incompressible elasticity equations deteriorates as the Lamé coefficient $λ\to\infty$, or equivalently as the Poisson ratio $ν\to1/2$. This phenomenon, known as locking or non-robustness, remains not fully understood despite extensive investigation. In this work, we illustrate first that an analogous instability arises when applying the popular Physics-Informed Neural Networks (PINNs) to nearly incompressible elasticity problems, leading to significant loss of accuracy and convergence difficulties. Then, to overcome this challenge, we propose a robust decomposition-based PINN framework that reformulates the elasticity equations into balanced subsystems, thereby eliminating the ill-conditioning that causes locking. Our approach simultaneously solves the forward and inverse problems to recover both the decomposed field variables and the associated external conditions. We will also perform a convergence analysis to further enhance the reliability of the proposed approach. Moreover, through various numerical experiments, including constant, variable and parametric Lamé coefficients, we illustrate the efficiency of the proposed methodology.
MLJan 31, 2020
Deep Learning Based Unsupervised and Semi-supervised Classification for KeratoconusNicole Hallett, Kai Yi, Josef Dick et al.
The transparent cornea is the window of the eye, facilitating the entry of light rays and controlling focusing the movement of the light within the eye. The cornea is critical, contributing to 75% of the refractive power of the eye. Keratoconus is a progressive and multifactorial corneal degenerative disease affecting 1 in 2000 individuals worldwide. Currently, there is no cure for keratoconus other than corneal transplantation for advanced stage keratoconus or corneal cross-linking, which can only halt KC progression. The ability to accurately identify subtle KC or KC progression is of vital clinical significance. To date, there has been little consensus on a useful model to classify KC patients, which therefore inhibits the ability to predict disease progression accurately. In this paper, we utilised machine learning to analyse data from 124 KC patients, including topographical and clinical variables. Both supervised multilayer perceptron and unsupervised variational autoencoder models were used to classify KC patients with reference to the existing Amsler-Krumeich (A-K) classification system. Both methods result in high accuracy, with the unsupervised method showing better performance. The result showed that the unsupervised method with a selection of 29 variables could be a powerful tool to provide an automatic classification tool for clinicians. These outcomes provide a platform for additional analysis for the progression and treatment of keratoconus.
NAJul 26, 2017
Covering of spheres by spherical caps and worst-case error for equal weight cubature in Sobolev spacesJohann S. Brauchart, Josef Dick, Edward B. Saff et al.
We prove that the covering radius of an $N$-point subset $X_N$ of the unit sphere $S^d \subset R^{d+1}$ is bounded above by a power of the worst-case error for equal weight cubature $\frac{1}{N}\sum_{\mathbf{x} \in X_N}f(\mathbf{x}) \approx \int_{S^d} f \, \mathrm{d} σ_d$ for functions in the Sobolev space $\mathbb{W}_p^s(S^d)$, where $σ_d$ denotes normalized area measure on $S^d.$ These bounds are close to optimal when $s$ is close to $d/p$. Our study of the worst-case error along with results of Brandolini et al. motivate the definition of Quasi-Monte Carlo (QMC) design sequences for $\mathbb{W}_p^s(S^d)$, which have previously been introduced only in the Hilbert space setting $p=2$. We say that a sequence $(X_N)$ of $N$-point configurations is a QMC-design sequence for $\mathbb{W}_p^s(S^d)$ with $s > d/p$ provided the worst-case equal weight cubature error for $X_N$ has order $N^{-s/d}$ as $N \to \infty$, a property that holds, in particular, for a sequence of spherical $t$-designs in which each design has order $t^d$ points. For the case $p = 1$, we deduce that any QMC-design sequence $(X_N)$ for $\mathbb{W}_1^s(S^d)$ with $s > d$ has the optimal covering property; i.e., the covering radius of $X_N$ has order $N^{-1/d}$ as $N \to \infty$. A significant portion of our effort is devoted to the formulation of the worst-case error in terms of a Bessel kernel, and showing that this kernel satisfies a Bernstein type inequality involving the mesh ratio of $X_N$. As a consequence we prove that any QMC-design sequence for $\mathbb{W}_p^s(S^d)$ is also a QMC-design sequence for $\mathbb{W}_{p^\prime}^s(S^d)$ for all $1 \leq p < p^\prime \leq \infty$ and, furthermore, if $(X_N)$ is a quasi-uniform QMC-design sequence for $\mathbb{W}_p^s(S^d)$, then it is also a QMC-design sequence for $\mathbb{W}_p^{s^\prime}(S^d)$ for all $s > s^\prime > d/p$.
NAAug 30, 2016
Lattice based integration algorithms: Kronecker sequences and rank-1 latticesJosef Dick, Friedrich Pillichshammer, Kosuke Suzuki et al.
We prove upper bounds on the order of convergence of lattice based algorithms 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 spaces of dominating mixed smoothness $\mathring{\mathbf{B}}^s_{p,θ}$, which also comprise the concept of Sobolev spaces of dominating mixed smoothness $\mathring{\mathbf{H}}^s_{p}$ as special cases. The considered algorithms are quasi-Monte Carlo rules with underlying nodes from $T_N(\mathbb{Z}^d) \cap [0,1)^d$, where $T_N$ is a real invertible generator matrix of size $d$. For such rules the worst-case error can be bounded in terms of the Zaremba index of the lattice $\mathbb{X}_N=T_N(\mathbb{Z}^d)$. We apply this result to Kronecker lattices and to rank-1 lattice point sets, which both lead to optimal error bounds up to $\log N$-factors for arbitrary smoothness $s$. The advantage of Kronecker lattices and classical lattice point sets is that the run-time of algorithms generating these point sets is very short.
NAAug 9, 2015
Multi-level higher order QMC Galerkin discretization for affine parametric operator equationsJosef Dick, Frances Kuo, Quoc Thong Le Gia et al.
We develop a convergence analysis of a multi-level algorithm combining higher order quasi-Monte Carlo (QMC) quadratures with general Petrov-Galerkin discretizations of countably affine parametric operator equations of elliptic and parabolic type, extending both the multi-level first order analysis in [\emph{F.Y.~Kuo, Ch.~Schwab, and I.H.~Sloan, Multi-level quasi-Monte Carlo finite element methods for a class of elliptic partial differential equations with random coefficient} (in review)] and the single level higher order analysis in [\emph{J.~Dick, F.Y.~Kuo, Q.T.~Le~Gia, D.~Nuyens, and Ch.~Schwab, Higher order QMC Galerkin discretization for parametric operator equations} (in review)]. We cover, in particular, both definite as well as indefinite, strongly elliptic systems of partial differential equations (PDEs) in non-smooth domains, and discuss in detail the impact of higher order derivatives of {\KL} eigenfunctions in the parametrization of random PDE inputs on the convergence results. Based on our \emph{a-priori} error bounds, concrete choices of algorithm parameters are proposed in order to achieve a prescribed accuracy under minimal computational work. Problem classes and sufficient conditions on data are identified where multi-level higher order QMC Petrov-Galerkin algorithms outperform the corresponding single level versions of these algorithms. Numerical experiments confirm the theoretical results.
NAJun 26, 2015
On a projection-corrected component-by-component constructionJosef Dick, Peter Kritzer
The component-by-component construction is the standard method of finding good lattice rules or polynomial lattice rules for numerical integration. Several authors have reported that in numerical experiments the generating vector sometimes has repeated components. We study a variation of the classical component-by-component algorithm for the construction of lattice or polynomial lattice point sets where the components are forced to differ from each other. This avoids the problem of having projections where all quadrature points lie on the main diagonal. Since the previous results on the worst-case error do not apply to this modified algorithm, we prove such an error bound here. We also discuss further restrictions on the choice of components in the component-by-component algorithm.
NAJun 24, 2015
Higher order Quasi-Monte Carlo integration for holomorphic, parametric operator equationsJosef Dick, Quoc T. Le Gia, Christoph Schwab
We analyze the convergence of higher order Quasi-Monte Carlo (QMC) quadratures of solution-functionals to countably-parametric, nonlinear operator equations with distributed uncertain parameters taking values in a separable Banach space $X$ admitting an unconditional Schauder basis. Such equations arise in numerical uncertainty quantification with random field inputs. Unconditional bases of $X$ render the random inputs and the solutions of the forward problem countably parametric, deterministic. We show that these parametric solutions belong to a class of weighted Bochner spaces of functions of countably many variables, with a particular structure of the QMC quadrature weights: up to a (problem-dependent, and possibly large) finite dimension, product weights can be used, and beyond this dimension, weighted spaces with so-called SPOD weights recently introduced in [F.Y.~Kuo, Ch.~Schwab, I.H.~Sloan, Quasi-Monte Carlo finite element methods for a class of elliptic partial differential equations with random coefficients. SIAM J. Numer. Anal., 50, 3351--3374, 2012.] can be used to describe the solution regularity. The regularity results in the present paper extend those in [J. Dick, F.Y.~Kuo, Q.T.~Le Gia, D.~Nuyens, Ch.~Schwab, Higher order QMC (Petrov-)Galerkin discretization for parametric operator equations. SIAM J. Numer. Anal., 52, 2676 -- 2702, 2014.] established for affine parametric, linear operator families; they imply, in particular, efficient constructions of (sequences of) QMC quadrature methods there, which are applicable to these problem classes. We present a hybridized version of the fast component-by-component (CBC for short) construction of a certain type of higher order digital net.
CODec 2, 2014
Discrepancy estimates for variance bounding Markov chain quasi-Monte CarloJosef Dick, Daniel Rudolf
Markov chain Monte Carlo (MCMC) simulations are modeled as driven by true random numbers. We consider variance bounding Markov chains driven by a deterministic sequence of numbers. The star-discrepancy provides a measure of efficiency of such Markov chain quasi-Monte Carlo methods. We define a pull-back discrepancy of the driver sequence and state a close relation to the star-discrepancy of the Markov chain-quasi Monte Carlo samples. We prove that there exists a deterministic driver sequence such that the discrepancies decrease almost with the Monte Carlo rate $n^{1/2}$. As for MCMC simulations, a burn-in period can also be taken into account for Markov chain quasi-Monte Carlo to reduce the influence of the initial state. In particular, our discrepancy bound leads to an estimate of the error for the computation of expectations. To illustrate our theory we provide an example for the Metropolis algorithm based on a ball walk. Furthermore, under additional assumptions we prove the existence of a driver sequence such that the discrepancy of the corresponding deterministic Markov chain sample decreases with order $n^{-1+δ}$ for every $δ>0$.
NANov 11, 2014
Numerical integration in $\log$-Korobov and $\log$-cosine spacesJosef Dick, Peter Kritzer, Gunther Leobacher et al.
QMC rules are equal weight quadrature rules for approximating integrals over $[0,1]^s$. One line of research studies the integration error of functions in the unit ball of so-called Korobov spaces, which are Hilbert spaces of periodic functions on $[0,1]^s$ with square integrable partial mixed derivatives of order $α$. Using Parseval's identity, this smoothness can be defined for all real numbers $α> 1/2$. This condition is necessary as otherwise the Korobov space contains discontinuous functions for which function evaluation is not well defined. This paper is concerned with more precise endpoint estimates of the integration error using QMC rules for Korobov spaces with $α$ arbitrarily close to $1/2$. To obtain such estimates we introduce a $\log$-scale for functions with smoothness close to $1/2$, which we call $\log$-Korobov spaces. We show that lattice rules can be used to obtain an integration error of order $\mathcal{O}(N^{-1/2} (\log N)^{-μ(1-λ)/2})$ for any $1/μ<λ\le 1$, where $μ>1$ is a power in the $\log$-scale. We also consider tractability of numerical integration for weighted Korobov spaces with product weights $(γ_j)_{j \in \mathbb{N}}$. It is known that if $\sum_{j=1}^\infty γ_j^τ< \infty$ for some $1/(2α) < τ\le 1$ one can obtain error bounds which are independent of the dimension. In this paper we give a more refined estimate for the case where $τ$ is close to $1/(2 α)$, namely we show dimension independent error bounds under the condition that $\sum_{j=1}^\infty γ_j \max\{1, \log γ_j^{-1}\}^{μ(1-λ)} < \infty$ for some $1/μ< λ\le 1$. The essential tool in our analysis is a $\log$-scale Jensen's inequality. The results described above also apply to integration in $\log$-cosine spaces using tent-transformed lattice rules.
CAOct 8, 2014
Functions of bounded variation, signed measures, and a general Koksma-Hlawka inequalityChristoph Aistleitner, Josef Dick
In this paper we prove a correspondence principle between multivariate functions of bounded variation in the sense of Hardy and Krause and signed measures of finite total variation, which allows us to obtain a simple proof of a generalized Koksma--Hlawka inequality for non-uniform measures. Applications of this inequality to importance sampling in Quasi-Monte Carlo integration and tractability theory are given. Furthermore, we discuss the problem of transforming a low-discrepancy sequence with respect to the uniform measure into a sequence with low discrepancy with respect to a general measure $μ$, and show the limitations of a method suggested by Chelson.