Adaptive Sparse Grids

Sparse grids can be adaptively constructed using a Gerstner–Griebel [4] type construction. This is implemented in the adaptive_sparsegrid function.

One dimensional example

To illustrate the adaptive algorithm consider functions that can be exactly approximated using a sparse grid approximation. Consider the polynomial $f:[-1,1]\mapsto\R$

\[f(x) = 3x^2 + 2x +1.\]

using SparseGridsKit
nparams = 1
f(x) = @. 3.0*x[1]^2 + 2.0*x[1] + 1.0
f (generic function with 1 method)

The adaptive_sparsegrid function can be called with a function $f$ and the dimension of the function domain nparams.

(sg, f_on_Z) = adaptive_sparsegrid(f, nparams)
SparseGridApproximation(sg, f_on_Z)
SparseGridApproximation(SparseGrid(1, [[-1.0, 1.0]], [1 2], [0, 1], Real[1.0; 0.0; -1.0;;], [0.16666666666666669, 0.6666666666666667, 0.16666666666666669], Points[CCPoints{Float64}([-1.0, 1.0])], Level[Doubling()], SparseGridsKit.sparse_grid_data{1}(SparseGridsKit.SparseTerm{1}[SparseGridsKit.SparseTerm{1}(([2, 1],)), SparseGridsKit.SparseTerm{1}(([2, 2],)), SparseGridsKit.SparseTerm{1}(([2, 3],))], [(2,), (1,), (3,)], [1, 2, 3], [1, 1, 1], [Dict([2, 3] => 3, [2, 1] => 2, [1, 1] => 1, [2, 2] => 1)], [[[0.0], [1.0, 0.0, -1.0]]], [[0.0, 1.0, -1.0]], [[[1.0], [0.16666666666666669, 0.6666666666666667, 0.16666666666666669]]])), [6.0, 1.0, 2.0], [[-1.0, 1.0]])

The function f can be exactly represented by a three point interpolant, which is identified in the adaptive algorithm.

using LinearAlgebra
test_points = range(-1,stop=1,length=100)
error = norm([[f(x)] - interpolate_on_sparsegrid(sg,f_on_Z,x) for x in test_points])
println("Error: $error")
Error: 3.6502539264068535e-15
get_n_grid_points(sg)
3

Taking powers $k$ of the polynomial $f$ gives a polynomial $f^k$ of polynomial degree $2^k$. Consequently, this is represented exactly using $2^k+1$ interpolation points. The adaptive algorithm identifies this.

f2(x) = f(x).^2
(sg, f2_on_Z) = adaptive_sparsegrid(f2, nparams)
# Expect three point approx cubic (1 iteration to suffice)
error = norm([f2(x)] ≈ interpolate_on_sparsegrid(sg,f2_on_Z,x) for x in test_points),
println("Error: $error")
(10.0, nothing)
get_n_grid_points(sg) == 5
true
f3(x) = f(x).^3
(sg, f3_on_Z) = adaptive_sparsegrid(f3, nparams)
# Expect three point approx cubic (1 iteration to suffice)
error = norm([f3(x)] ≈ interpolate_on_sparsegrid(sg,f3_on_Z,x) for x in test_points)
println("Error: $error")
Error: 10.0
get_n_grid_points(sg)
9

Multi-dimensional example

The adaptive algorithm is also applicable for higher dimensional functions. Consider the genz gaussianpeak example.

using SparseGridsKit, LinearAlgebra
n = 8
C = 1.0
W = 0.0
T = "exponentialdecay"
N = "gaussianpeak"
f = genz(n::Int, C::Float64, W::Float64, T::String, N::String)
#5 (generic function with 1 method)

An adaptive approximation is formed. This produces a (relatively) accurate approximation of the high-dimensional function $f$. The iterations halt when the stopping criterion of max profit p_max < proftol.

# Approximate
(sg, f_on_Z) = adaptive_sparsegrid(f, n)

nmc = 100
test_points = [2 * rand(8) .- 1 for ii = 1:nmc]

f_approx_on_test = interpolate_on_sparsegrid(sg,f_on_Z,test_points)

norm([f(x) - f_approx_on_test[i] for (i,x) in enumerate(test_points)])
0.0015471677772621146

Reducing proftol allows further iterations and a more accurate approximation.

(sg, f_on_Z) = adaptive_sparsegrid(f, n, proftol=1e-5)

f_approx_on_test = interpolate_on_sparsegrid(sg,f_on_Z,test_points)

norm([f(x) - f_approx_on_test[i] for (i,x) in enumerate(test_points)])
0.00041235982884738363

Function Reference

SparseGridsKit.EvaluationDictionaryMethod
EvaluationDictionary(sg,f_on_z)

Creates an EvaluationDictionary from a sparse grid and function evaluations

Arguments

  • sg: Sparse grid
  • f_on_z: Function evaluations on the sparse grid

Returns

  • evaluations: EvaluationDictionary
source
SparseGridsKit.adaptive_estimateMethod
adaptive_estimate(sg, datastore, pcl, type=:deltaint)

Estimates the profit of adding multi-indices {α} in reduced margin to the sparse grid.

Arguments

  • sg: Sparse grid.
  • datastore: EvaluationDictionary.
  • rule: Level function(s) for sparse grid.
  • knots: Knot function(s) for sparse grid.
  • type: Type of profit computation. Default is :deltaint.

Returns

  • Vector of profits for each multi-index α.
source
SparseGridsKit.adaptive_markFunction
adaptive_mark(α, p_α, θ=1e-4)

Dorfler style marking of multi-indices based on profits

Arguments

  • α: Multi-indices in reduced margin.
  • p_α: Vector of profits.
  • θ: Threshold for marking. Default is 1e-4.
source
SparseGridsKit.adaptive_refineMethod
adaptive_refine(sg, α_marked, rule, knots)

Refines the sparse grid based on marked multi-indices

Arguments

  • sg: Sparse grid.
  • α_marked: Marked multi-indices.
  • rule: Level function(s) for sparse grid.
  • knots: Knot function(s) for sparse grid.

Returns

  • sg: Refined sparse grid.
  • f_on_z: Function evaluations on the refined sparse grid.
source
SparseGridsKit.adaptive_solve!Method
adaptive_solve!(f, datastore, sg)

Adds function evaluations on sg to datastore

Arguments

  • f: Function
  • datastore: EvaluationDictionary
  • sg: Sparse grid

Returns

  • datastore: Updated EvaluationDictionary
source
SparseGridsKit.adaptive_sparsegridMethod
adaptive_sparsegrid(f, nparams; maxpts = 100, proftol = 1e-4, rule=Doubling(), knots=CCPoints(), θ=1e-4, type=:deltaint)

Constructs an adaptive sparse grid for approximating the function f in nparams dimensions.

Arguments

  • f: Function to be approximated taking arguments x with length(x)=nparams.
  • nparams: Dimension of domain.
  • maxpts: (Optional) Maximum number of points to include in the sparse grid. Default is 100.
  • proftol: (Optional) Tolerance for profits. Default is 1e-4.
  • rule: (Optional) Level function(s) for sparse grid. Default is Doubling().
  • knots: (Optional) Knot function(s) for sparse grid. Default is CCPoints().
  • θ: (Optional) Threshold for marking. Default is 1e-4.
  • type: (Optional) Type of profit computation. Default is :deltaint.
  • costfunction (Optional) Cost function for fidelity related multi-indices.

Returns

  • sg: Final sparse grid used to approximate f
  • f_on_z: Evaluations of f on grid points in sparse grid sg
source
SparseGridsKit.adaptive_terminateMethod

adaptiveterminate(sg, pα, maxpts, proftol)

Test profits and sparse grid to determine loop termination

Arguments

  • sg: Sparse grid.
  • p_α: Vector of profits.
  • maxpts: Maximum number of grid points allowed.
  • proftol: Profit tolerance.

Returns

  • Flag to indicate loop termination if maxpts has been reached, or proftol attained.
source
SparseGridsKit.add_evaluation!Method
add_evaluation!(evaluations, z, f_at_z)

Adds the evaluation fatz into the evaluations_database

Arguments

  • evaluations: EvaluationDictionary
  • z: Point
  • f_at_z: Function evaluation at z

Returns

  • evaluations: Updated EvaluationDictionary
source
SparseGridsKit.add_evaluations!Method
add_evaluations!(evaluations, sg::SparseGrid, f_on_z)

Adds the evaluations fonz into the evaluations_database

Arguments

  • evaluations: EvaluationDictionary
  • sg: Sparse grid
  • f_on_z: Function evaluations on the sparse grid

Returns

  • evaluations: Updated EvaluationDictionary
source
SparseGridsKit.compute_profitMethod
compute_profit(sg_α, f_α, f, cost, pcl; type=:deltaint)

Computes the "profit" of a sparse grid supplemented by a multi-index α

Arguments

  • sg_α: Enhanced sparse grid with α
  • sg: Sparse grid
  • f_α: Function evaluations at grid points in α
  • f: Function evaluations on the sparse grid.
  • cost: Cost of adding α to the sparse grid.
  • type: Type of profit computation. Default is :deltaint. Options are :deltaint, :deltaintcost, :Linf, :Linfcost.

Returns

  • Computed profit as Expected change in approximation.
source
SparseGridsKit.retrieve_evaluationsMethod
retrieve_evaluations(evaluations, sg)

Retrieves the evaluations from the evaluations_database

Arguments

  • evaluations: EvaluationDictionary
  • sg: Sparse grid

Returns

  • f_on_z: Function evaluations on the sparse grid
source