How do error estimates scale for multidimensional cubature?

Question

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.

Answer 1

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})$$