SPARSE_GRID_HW
Sparse Grids for Uniform and Normal Weights
Heiss and Winschel

SPARSE_GRID_HW is a FORTRAN77 library which can compute sparse grids for multidimensional integration, based on 1D rules for the unit interval with unit weight function, or for the real line with the Gauss-Hermite weight function. The original MATLAB code is by Florian Heiss and Viktor Winschel.

The FORTRAN77 version of the library is only partially complete. Completion of the translation has been held up by the fact that FORTRAN77 does not "officially" support allocatable arrays, and so the sizes of certain arrays must be determined and allocated in advance.

The original version of this software, and other information, is available at http://sparse-grids.de .

Four built-in 1D families of quadrature rules are supplied, and the user can extend the package by supplying any family of 1D quadrature rules.

The built-in families are identified by a 3-letter key which is also the name of the FORTRAN routine that returns members of the family:

• gqu, standard Gauss-Legendre quadrature rules, for the unit interval [0,1], with weight function w(x) = 1.
• gqn, standard Gauss-Hermite quadrature rules, for the infinite interval (-oo,+oo), with weight function w(x) = exp(-x*x/2)/sqrt(2*pi).
• kpu, Kronrod-Patterson quadrature rules, for the unit interval [0,1], with weight function w(x) = 1. These sacrifice some of the precision of gqu in order to provide a family of nested rules.
• kpn, Kronrod-Patterson quadrature rules, for the infinite interval (-oo,+oo), with weight function w(x) = exp(-x*x/2)/sqrt(2*pi). These sacrifice some of the precision of gqn in order to provide a family of nested rules.

The user can build new sparse grids by supplying a 1D quadrature family. Examples provided include:

• cce_order, Clenshaw-Curtis Exponential quadrature rules, for the unit interval [0,1], with weight function w(x) = 1. The K-th call returns the rule of order 1 if K is 1, and 2*(K-1)+1 otherwise.
• ccl_order, Clenshaw-Curtis Linear quadrature rules, for the unit interval [0,1], with weight function w(x) = 1. The K-th call returns the rule of order 2*K-1.
• ccs_order, slow Clenshaw-Curtis Slow quadrature rules, for the unit interval [0,1], with weight function w(x) = 1. The K-th call returns the rule of order 1 if K is 1, and otherwise a rule whose order N has the form 2^E+1 and is the lowest such order with precision at least 2*K-1.

Licensing:

The computer code and data files described and made available on this web page are distributed under the GNU LGPL license.

Languages:

SPARSE_GRID_HW is available in a C version and a C++ version and a FORTRAN77 version and a FORTRAN90 version and a MATLAB version

Related Data and Programs:

GRID_DISPLAY, a MATLAB library which can display a 2D or 3D grid or sparse grid.

NINT_EXACTNESS_MIXED, a FORTRAN90 program which measures the polynomial exactness of a multidimensional quadrature rule based on a mixture of 1D quadrature rule factors.

PRODUCT_RULE, a FORTRAN90 program which constructs a product quadrature rule from identical 1D factor rules.

QUADRULE, a FORTRAN77 library which defines quadrature rules for various intervals and weight functions.

SANDIA_RULES, a FORTRAN90 library which generates Gauss quadrature rules of various orders and types.

SANDIA_SPARSE, a FORTRAN90 library which computes the points and weights of a Smolyak sparse grid, based on a variety of 1-dimensional quadrature rules.

SMOLPACK, a C library which implements Novak and Ritter's method for estimating the integral of a function over a multidimensional hypercube using sparse grids, by Knut Petras.

SPARSE_GRID_HERMITE, a FORTRAN90 library which creates sparse grids based on Gauss-Hermite rules.

SPINTERP, a MATLAB library which carries out piecewise multilinear hierarchical sparse grid interpolation; an earlier version of this software is ACM TOMS Algorithm 847, by Andreas Klimke;

Author:

Original MATLAB code by Florian Heiss and Viktor Winschel. FORTRAN77 version by John Burkardt.

Reference:

• Alan Genz, Bradley Keister,
  Fully symmetric interpolatory rules for multiple integrals over infinite regions with Gaussian weight,
  Journal of Computational and Applied Mathematics,
  Volume 71, 1996, pages 299-309.
• Florian Heiss, Viktor Winschel,
  Likelihood approximation by numerical integration on sparse grids,
  Journal of Econometrics,
  Volume 144, Number 1, May 2008, pages 62-80.
• Thomas Patterson,
  The optimal addition of points to quadrature formulae,
  Mathematics of Computation,
  Volume 22, Number 104, October 1968, pages 847-856.
• Knut Petras,
  Smolyak Cubature of Given Polynomial Degree with Few Nodes for Increasing Dimension,
  Numerische Mathematik,
  Volume 93, Number 4, February 2003, pages 729-753.

Source Code:

• sparse_grid_hw.f, the source code.
• sparse_grid_hw.sh, BASH commands to compile the source code.

Examples and Tests:

• sparse_grid_hw_prb.f, a sample calling program.
• sparse_grid_hw_prb.sh, BASH commands to compile and run the sample program.
• sparse_grid_hw_prb_output.txt, the output file.

List of Routines:

• CCE_ORDER: order of a Clenshaw Curtis Exponential rule from the level.
• CCL_ORDER: order of a Clenshaw Curtis Linear (CCL) quadrature rule.
• CCS_ORDER: order of a "slow growth" Clenshaw Curtis quadrature rule.
• CC computes a Clenshaw Curtis quadrature rule based on order.
• FN_INTEGRAL is the integral of the Hermite test function.
• FN_VALUE is a Hermite test function.
• FU_INTEGRAL is the integral of the test function for the [0,1]^D interval.
• FU_VALUE is a sample function for the [0,1]^D interval.
• GET_SEQ generates all positive integer D-vectors that sum to NORM.
• GQN provides data for Gauss quadrature with a normal weight.
• GQN_ORDER computes the order of a GQN rule from the level.
• GQN2_ORDER computes the order of a GQN rule from the level.
• GQU provides data for Gauss quadrature with a uniform weight.
• GQU_ORDER computes the order of a GQU rule from the level.
• I4_CHOOSE computes the binomial coefficient C(N,K).
• I4_FACTORIAL2 computes the double factorial function.
• I4_MOP returns the I-th power of -1 as an I4 value.
• I4MAT_PRINT prints an I4MAT.
• I4MAT_PRINT_SOME prints some of an I4MAT.
• I4VEC_CUM0 computes the cumulutive sum of the entries of an I4VEC.
• I4VEC_PRINT prints an I4VEC.
• I4VEC_PRODUCT returns the product of the entries of an I4VEC.
• I4VEC_SUM returns the sum of the entries of an I4VEC.
• I4VEC_TRANSPOSE_PRINT prints an I4VEC "transposed".
• KPN provides data for Kronrod-Patterson quadrature with a normal weight.
• KPN_ORDER computes the order of a KPN rule from the level.
• KPU provides data for Kronrod-Patterson quadrature with a uniform weight.
• KPU_ORDER computes the order of a KPU rule from the level.
• NUM_SEQ returns the number of compositions of the integer N into K parts.
• NWSPGR_SIZE determines the size of a sparse grid rule.
• QUAD_RULE_PRINT prints a multidimensional quadrature rule.
• R8_UNIFORM_01 returns a pseudorandom R8 scaled to [0,1].
• R8CVV_OFFSET determines the row offsets of an R8CVV.
• R8CVV_PRINT prints an R8CVV.
• R8CVV_RGET gets row I from an R8CVV.
• R8CVV_RSET sets row I from an R8CVV.
• R8MAT_NORMAL_01 returns a unit pseudonormal R8MAT.
• R8MAT_TRANSPOSE_PRINT prints an R8MAT, transposed.
• R8MAT_TRANSPOSE_PRINT_SOME prints some of an R8MAT transposed.
• R8MAT_UNIFORM_01 returns a unit pseudorandom R8MAT.
• R8VEC_COPY copies an R8VEC.
• R8VEC_DIRECT_PRODUCT creates a direct product of R8VEC's.
• R8VEC_DIRECT_PRODUCT2 creates a direct product of R8VEC's.
• R8VEC_DOT_PRODUCT finds the dot product of a pair of R8VEC's.
• R8VEC_NORMAL_01 returns a unit pseudonormal R8VEC.
• R8VEC_PRINT prints an R8VEC.
• R8VEC_SUM sums the entries of an R8VEC.
• R8VEC_UNIFORM_01 returns a unit pseudorandom R8VEC.
• RULE_ADJUST adjusts a 1D quadrature rule from [A,B] to [C,D].
• RULE_SORT sorts a multidimensional quadrature rule.
• SORT_HEAP_EXTERNAL externally sorts a list of items into ascending order.
• SYMMETRIC_SPARSE_SIZE sizes a symmetric sparse rule.
• TENSOR_PRODUCT generates a tensor product quadrature rule.
• TENSOR_PRODUCT_CELL generates a tensor product quadrature rule.
• TIMESTAMP prints out the current YMDHMS date as a timestamp.

You can go up one level to the FORTRAN77 source codes.