# How do error estimates scale for multidimensional cubature?

Suppose I have a quadrature method with a theoretical error estimate that scales with the number of points as $e(N)$:

$$\int_a^b f(x)\mathrm{d}x = \sum_{i=1}^N w_i f(x_i) + \mathcal{O}\bigl(e(N)\bigr)$$

For example, $e(N) = N^{-4}$ for Simpson's rule, or $e(N) = N^{-2}$ for the trapezoidal rule. In this question I only care about how the error estimate scales with $N$, not its dependence on properties of the function or any leading coefficients. (And yes, I know it's just an estimate.)

If I use this quadrature method to perform a $d$-dimensional integral by repeated 1D integration,

$$\int_{a_0}^{b_0} \int_{a_1(x_0)}^{b_1(x_0)} \cdots \int_{a_d(x_0,\ldots,x_{d-1})}^{b_d(x_0,\ldots,x_{d-1})} f(\vec{x})\mathrm{d}^d\vec{x}$$

how does the error estimate for the multidimensional integral scale with $N$, the number of function evaluations along each dimension? Is it just $e(N)$? Or some function of $e(N)$? (like $[e(N)]^d$) Is it impossible to give a general relation without knowing the specific quadrature method in use?

I've done some Google searching and checking in Numerical Recipes but I can't find a straightforward answer to this. I would have thought it would be common.

Using product quadrature rules for multi-dimensional integrals suffers from the so-called curse of dimensionality. An $O(N^{-2})$ accurate rule using N evaluations in one-dimension is generally $O(N^{-2})$ accurate when applied as a product rule in multiple dimensions, but there will be $M = N^d$ evaluations required. So the accuracy is $O(M^{-2/d}).$ The required number of evaluations becomes intractable beyond 10 dimensions.

Monte Carlo and Quasi-Monte Carlo methods seem to be preferred approaches for many dimensions with an error that diminishes like $1/\sqrt{M}$.

For example, if $f:[0,1] \rightarrow \mathbb{R}$ and the second-derivative is bounded -- $|f''(x)|\leq K$ -- then the error bound for the mid-point rule applied at the mid-points $\bar{x_i}= (2i-1)/(2N)$ of $N$ evenly spaced points $x_i = i/N$ is

$$\left|\int_{0}^{1}f(x)\,dx - \frac1{N}\sum_{i=1}^{N}f\left(\bar{x_i}\right)\right|\leq \frac{K}{24N^2}.$$

In the two-dimensional case, $f:[0,1]^2 \rightarrow \mathbb{R}$, where $f$ has bounded second-order partial derivatives we have the same order of accuracy with respect to N. Using a Taylor series approximation, we find that the first- and mixed second-order derivative terms vanish after integrating and the remaining error is

$$\left|\int_{0}^{1}\int_{0}^{1}f(x,y)\,dx\,dy - \frac1{N^2}\sum_{i=1}^{N}\sum_{j=1}^{N}f\left(\bar{x_i},\bar{y_i}\right)\right|\\=\left|\sum_{i=1}^{N}\sum_{j=1}^{N}\int_{x_{i-1}}^{x_i}\int_{y_{j-1}}^{y_j}\frac1{2}\left[f_{xx}(\xi_x,\xi_y)(x-\bar{x_i})^2+f_{yy}(\eta_x,\eta_y)(y-\bar{y_i})^2\right]\,dx\,dy \right| \\ = N^2O(N^{-4})= O(N^{-2})$$

• Makes sense. And the same argument should hold for algorithms like Gaussian quadrature with exponential convergence $e(N)\sim a^N$, right? Jul 16, 2014 at 21:40
• @Davd Z. Yes -- I did a quick check and that convergence rate carries over to more dimensions. Intuitively extending to higher dimensions shouldn't shouldn't make it faster. Despite that improved rate competing against the exponential increase in evaluations with d it stills seems beyond 15 dimensions the cubature is not effective for most functions. Jul 17, 2014 at 1:23