Mathematical Background¶
This page is still a work in progress. I will try to provide a very basic description of each of the pieces that are necessary.
Smolyak Grid¶
One of the standard forms of building an approximation grid is to use a simple tensor-product. While this seems relatively trivial for low dimensions, but the number of points to evaluate quickly becomes intractable for larger dimensions. In Bellman (1961) this problem is referred to as the “curse of dimensionality.” Two years after Bellman’s paper, Sergey Smolyak introduced a numerical technique where the number of grid points needed to approximate grew polynomially instead of exponentially. As stated in Judd, Maliar, Maliar, Valero; the idea behind this technique was “that some elements produced by tensor-product rules are more important for representing multidimensional functions than the others.” The tensor-product typically takes as parameters the dimension of the grid and the number of points, \(n\) to be evaluated at in each dimension which produces a grid with \(n^d\) points (Note: There are other ways of doing this where the number of points is different in each dimension, but the resulting number of points is similar). The Smolyak grid takes an “accuracy” parameter \(\mu\) and the number of dimensions as parameters. The number of points in the grid is determined by the dimension and \(\mu\). The number of Smolyak grid points at \(\mu=1\) is \(1 + 2d\), at \(\mu=2\) it is \(1+4d+4d(d-1)\), etc... Notice that the number of grid points grows linearly at \(\mu=1\), and quadratically at \(\mu=2\).
The standard construction of a Smolyak grid uses nested sets of points. One typically uses the extrema of the Chebyshev Polynomials, which are known as the Chebyshev-Gauss-Lobatto points. We will continue our description using these points since the code is implemented using them. The nested sets require two conditions. First, that each set \(S_i\) has \(2^{i-1}\) points for \(i \geq 2\) and 1 if \(i=1\). Secondly, that the sets are nested. The first four nested sets are:
One then takes the tensor product of the unidimensional sets and then picks out the products that satisfy \(d \leq \sum_{j=1}^d i_j \leq d + \mu\) where \(i_j\) is the index of the unidimensional sets. For example if the parameters were \(d = 2, \mu = 2\) then you would build the first four nested sets of points (as shown above) and then take all of the tensor products that satisfied \(2 \leq i_1 + i_2 \leq 4\) which would give you \(S_1 \times S_1, S_1 \times S_2, S_2 \times S_1, S_2 \times S_2\). Then these would be your points. We will define the set of grid points to be \(\mathcal{H}^{d, \mu} := \bigcup_{d \leq |i^*| \leq d + \mu} \prod{i=1}^{d+\mu} S_i\) where \(i^* = \begin{bmatrix} i_1 & \dots i_d \end{bmatrix}\).
You can see how there would be repeated points and hence this method could be improved upon. This is one of the key results of the Judd, Maliar, Maliar, Valero (2013) paper. Instead of building nested sets, they build disjoint sets \(A_i\) such that \(A_1 = \{ 0 \}\) and \(A_i = S_i \backslash S_{i-1}\) for all \(i \geq 2\). Then the points are taken from the tensor-products of these sets in the same fashion as described above.
The following is a computationally efficient way of finding these grid points. The first step is to create the first \(\mu + 1\) unidimensional disjoint sets \(A_i\) as described above. Then one should create a vector of possible values i.e. \(\begin{bmatrix} 1, 2, \dots, \mu + 1 \end{bmatrix}\). It is important to note that these possible values only range from 1 to \(\mu + 1\) since the smallest possible index is 1 (To see this think of given \(d\) indexes. The smallest the first \(d-1\) of them can be is 1 which would sum to \(d-1\) hence the last index could only be valued up to \(\mu + 1\). We can then check all of the combinations (with replacement) of these to find the sets of indexes that would work. Once we have these we can permute them to capture all of the indexes that would work. Then you only take the tensor-products of these sets (Add reference to our code smol_inds and build_grid here).