sgpykit.main

Functions

`adapt_sparse_grid`(f, N_full, knots, ...)	Compute an adaptive sparse grid approximation of a function f.
`compute_sobol_indices_from_sparse_grid`(S, ...)	Compute Sobol indices from a sparse grid approximation of a scalar-valued function.
`convert_to_modal`(S, Sr, nodal_values, ...)	Convert a sparse grid interpolant to a modal (spectral) expansion.
`create_sparse_grid`(N, w, knots, lev2knots[, ...])	Generate a sparse grid (and corresponding quadrature weights) as a linear combination of full tensor grids, employing formula (2.9) of [Nobile-Tempone-Webster, SINUM 46/5, pages 2309-2345].
`create_sparse_grid_add_multiidx`(new_idx, ...)	Produce a grid obtained by adding a single multi-idx to a previously existing grid.
`create_sparse_grid_multiidx_set`(C, knots, ...)	Create a sparse grid starting from a multiindex-set rather than from a rule IDXSET(I) <= W.
`create_sparse_grid_quick_preset`(N, w)	Create a sparse grid with a quick preset configuration.
`derive_sparse_grid`(S, Sr, values_on_grid, ...)	Compute derivatives (gradients) of a scalar-valued function f: R^N -> R by centered finite differences formulas applied to the sparse grid approximation of f.
`evaluate_on_sparse_grid`(f, S, Sr[, ...])	Evaluate a function on a sparse grid, possibly recycling previous calls.
`export_sparse_grid_to_file`(Sr[, filename, ...])	Save knots of a reduced sparse grid to an ASCII file.
`hessian_sparse_grid`(S, Sr, values_on_grid, ...)	Compute the Hessian of a scalar-valued function f: R^N -> R by finite differences applied to the sparse grids approximation of f.
`interpolate_on_sparse_grid`(S, Sr, ...)	Interpolate a function on a sparse grid.
`plot3_sparse_grid`(S, dims, args, *kwargs)	Plot a sparse grid in 3D.
`plot_sparse_grid`(ax, S, dims, args, *kwargs)	Plot object based on the provided structure S and dimensions.
`plot_sparse_grids_interpolant`(ax, S, Sr, ...)	Plot the sparse grid interpolant of a function.
`quadrature_on_sparse_grid`(f, S, Sr[, ...])	Compute the integral of a function using a sparse grid.
`reduce_sparse_grid`(S[, tol])	Reduce a sparse grid by removing duplicate points within a given tolerance.

sgpykit.main.adapt_sparse_grid(f, N_full, knots, lev2knots, prev_adapt, controls)[source]

Compute an adaptive sparse grid approximation of a function f.

This function computes a sparse grid adapted to the approximation of a function F: Gamma subset R^N_full -> R^V according to some profit indicator. The implementation closely resembles the algorithm proposed in

Gerstner, M. Griebel, Dimension-Adaptive Tensor-Product Quadrature, Computing 71, 65-87 (2003)

and the various definitions of profits are discussed in

F. Nobile, L. Tamellini, F. Tesei, and R. Tempone, An Adaptive Sparse Grid Algorithm for Elliptic PDEs with Lognormal Diffusion Coefficient, Sparse Grids and Applications - Stuttgart 2014, Lecture Notes in Computational Science and Engineering 109, p. 191-220, Springer International Publishing Switzerland 2016

Roughly speaking, at each iteration the algorithm considers the multiidx with the highest profit among those already explored, adds it to the sparse grid, computes the profits of all its neighbouring indices (see controls.prof) and adds them to the list of explored indices.

Parameters:

fcallable

A function-handle, taking as input a column vector (N_full components) and returning a column-vector value (V components).

N_fullint

A scalar, denotes the dimension of the space Gamma over which F and the sparse grid are defined. Note that this is the full dimension, yet the algorithm might explore a smaller subset of variables to start with, N_curr, as specified in CONTROLS.var_buffer_size (see below)

knotscallable or list of callables

KNOTS can be a function handle or a cell arrays of function handles to differentiate between directions, but they must be either all nested or all non-nested formulae (the code cannot work partially with nested and partially with non-nested points) KNOTS cannot be used with the ‘nonprob’ option if later the “buffer” controls functionality is used, see later

lev2knotscallable

LEV2KNOTS must be a single function handle, i.e., it is not possible to use e.g., lev2knots_2step in one direction and lev2knots_lin in the next one.

For instance, in N=2, the following setups are ok:

knots = [lambda n: knots_uniform(n, 0, 1), lambda n: knots_uniform(n, 3, 5)], lev2knots = lev2knots_lin (same family on different intervals, same lev2knots)
knots = [lambda n: knots_CC(n, 0, 1), lambda n: knots_CC(n, 3, 5)], lev2knots = lev2knots_doubling (same family on different intervals, same lev2knots)

while the following ones are not:

knots = [lambda n: knots_uniform(n, 0, 1), lambda n: knots_CC(n, 3, 5)], lev2knots = lev2knots_doubling (one nested family, one non-nested family)
knots = [lambda n: knots_uniform(n, 0, 1), lambda n: knots_uniform(n, 3, 5)], lev2knots = [lev2knots_lin, lev2knots_doubling] (different lev2knots in different directions)
knots = [lambda n: knots_CC(n, 0, 1, 'nonprob'), lambda n: knots_CC(n, 3, 5, 'nonprob')], lev2knots = lev2knots_doubling (nonprob weights, cannot be used together with buffer)

See tutorial_adaptive.ipynb for an example that shows this.

prev_adaptdict or None

By setting PREV_ADAPT as the output ADAPTED of a previous computation, the new computation will resume from where the previous one left. Set PREV_ADAPT = None for a fresh start

controlsdict

A struct defining several parameters that control the algorithm flow. It is possible to define only some (or none) of these parameters, the others will be set to default value. Only the parameter ‘nested’ is mandatory. The parameter ‘pdf’ is mandatory only if controls.prof is set to any of the ‘weighted*’ choices (see below). For safety, ‘pdf’ is not set to default and an error will be raised if the field is not set controls.nested : bool

(mandatory) true if nested points will be used, false otherwise

controls.max_ptsint

is the maximum number of points the adapted sparse grid can have (default 1000). Observe that if the max nb of points is reached while checking the neighbours of the current idx, the algorithm will not be interrupted.

controls.prof_tolfloat

the algorithm will stop as soon as profit of the considered idx is smaller then prof_tol (default 1e-14)

controls.recyclingstr

can be set to ‘priority_to_evaluation’ (default) or ‘priority_to_recycling’. This affects only when non-nested points are used. Because in this case points enter and exit from the sparse grid, we should keep a record of all points that have been evaluated and recycle from it whenever possible. However, doing this exactly is very expensive in terms of computational cost with the current data structure if N is large. So, unless your function evaluation is very expensive, we recommend to leave the default setting, which does a “partial search”, i.e. searches for recycling only in the previous grid (instead than on the entire sparse grid history) and therefore implies that the same point might be evaluated multiple times.

controls.profstr

chooses how to compute the profit of an idx. In general, this entails comparing the sparse grids approximation of f before adding idx to the sparse grids (call this approximation S) and after having added it (call this approximation T). Note that S, T are actually functions of y in Gamma and return a vector with V components (same as f, the function to be approximated). In general, profit takes either of the two forms P(idx) = Op_y ( Op_vect(S,T) ) , P(idx) = Op_vect ( Op_y(S,T) ) where Op_vect acts as a “vector norm” over the V components of S,T (e.g., the euclidean norm for fixed y), and Op_y acts as a norm over Gamma (e.g. the maximum over a discrete set of values of y). Currently, the choices below are available. In all of these choices, Op_vect can be changed by setting the field controls.op_vect (see below)

‘Linf’ :: P(idx) = Op_y( Op_vect(S(y),T(y) ) where Op_vect(S(y),T(y)) = euclidean norm of S(y) - T(y) for a fixed y; Op_y = max over a discrete set of points of Gamma. For nested points, the algorithm considers the points that would be added to the grid by the considered idx, evalutes the function and the previous sparse grid on such points, and takes the max of such difference as profit estimate. For non-nested points, the sparse grid approx is not interpolant, hence we consider the difference between both the new and the previous sparse grids on the new points.
‘Linf/new_points’ (default) :: the ‘Linf’ profit estimate above is further divided by the number of new points considered (i.e. gain-cost ratio). For non-nested points, this is actually the cost of the tensor grid associated to the current midx (see Nobile Tamellini Tempone, ``Convergence of quasi optimal sparse grid approximation …’’
‘weighted Linf’ :: same as ‘Linf’, but the difference between previous sparse grid and function evaluation (or between new and old sparse grids, depending whether nested or non-nested points are used) is multiplied by the probability density function (see below for details). Useful when considering sparse grids on unbounded domains.
‘weighted Linf/new_points’ :: the ‘weighted Linf’ profit estimate is divided by the number of new points
‘deltaint’ :: P(idx) = Op_vect ( Op_y(S,T) ) = euclidean_norm( expct(S) - expct(T)) i.e. the profit estimate is the euclidean norm of the difference of the quadratures using the sparse grids with and without the given idx
‘deltaint/new_points’ :: the ‘deltaint’ profit estimate is divided by the number of new points

controls.op_vectcallable

changes the operator op_vect above. It is defined as a handle function lambda A, B: some_function_of_(A, B) where A B are two matrices with V rows, containing the evaluation of the two operators S T discussed above as columns, e.g. A=[S(y1) S(y2) S(y3) …], B=[T(y1) T(y2) T(y3) …]

controls.pts_tolfloat

is the tolerance used by evaluate_on_sparse_grid to check for equal points (default 1e-14)

controls.pdfcallable

pdf over Gamma subset R^N_full, the set over which the sparse of the random variables, to be used in weighted* profits. It must be provided as a function handle that takes as input a matrix with N rows (for generic N), i.e. where each column is a point where the pdf must be evaluated. For instance, if the weight is the standard gaussian pdf on each direction, controls.pdf = lambda Y: np.prod(norm.pdf(Y, 0, 1), axis=0), with Y matrix with N rows. For safety, ‘pdf’ is not set to default and an error will be raised if the field is not set

controls.var_buffer_sizeint

The algorithm starts exploring N_curr = <var_buffer_size> dimensions. As soon as points are placed in one dimension, a new dimension is added to the set of explored variables, i.e., N_curr = N_curr+1. In this way we ensure that there are always <var_buffer_size> explored but “non-activated” variables, i.e., along which no point is placed (default N_full).

IMPORTANT NOTE: Use this option with care! Is evaluating the function at a subset of variables equivalent to evaluating the function at the entire set of variables, where the ones outside the buffer are fixed at the mid-point of the interval? If this is not the case, the buffer cannot be used.

For instance:

f1 = 1/exp(y1*c1 + y2*c2 + y3*c3), with y1,y2 in [-0.5,0.5], y3 in [-0.2,0.2] This function is OK because f1 restricted to y1,y2 is y1*c1 + y2*c2 = y1*c1 + y2*c2 + 0*c3 = y1*c1 + y2*c2 + y3*c3 with y3 held at midpoint
f2 = 1/exp(y1*c1 + y2*c2 + y3*c3), with y1,y2,y3 in [0,1] (note the interval not centered in 0) This function is NOT OK because f1 restricted to y1,y2 is y1*c1 + y2*c2 ≠ y1*c1 + y2*c2 + 0.5*c3 = y1*c1 + y2*c2 + y3*c3 with y3 held at midpoint
f3 = y1 * y2 * y3, with y1,y2,y3 in [0,1] Is NOT OK because f1 restricted to y1,y2 is y1*y2 ≠ y1*y2*0.5 = y1*y2*y3 with y3 held at midpoint
f4 = cos(y1) + cos(y2) + cos(y3), with y1,y2,y3 in [-1,1] Is NOT OK because f1 restricted to y1,y2 is cos(y1) + cos(y2) ≠ cos(y1) + cos(y2) + 1 = cos(y1) + cos(y2) + cos(y3) with y3 held at midpoint
f5 = y1 * y2 * y3, with y1,y2,y3 in [0,2] Is OK because f1 restricted to y1,y2 is y1*y2 = y1*y2*1 = y1*y2*y3 with y3 held at midpoint

See tutorial_adaptive.ipynb for an example that shows this.

controls.plotbool

plot multiidx set and pause (default false)

Returns:

adapteddict

A struct containing the results of the algorithms. The same structure can be passed as input (PREV_ADAPT) to resume the computation. Its fields are adapted.S : the adapted sparse grid (built over both all indices whose profit has been computed,

even those whose profit hasn’t yet been chosen as the best one, i.e. those whose neighbours haven’t been explored yet).

adapted.Sr : its reduced version; adapted.f_on_Sr : the evaluations of f over the list of points contained in Sr; adapted.nb_pts : the number of points in Sr; adapted.nested : true if points used are nested adapted.nb_pts_visited : the number of points visited while building the sparse grid. For non-nested

points, this will be larger than nb_pts, because some points enter and then exit the grid when the corresponding idx exits from the combination technique.

adapted.num_evalsthe number of function evaluations done to compute the sparse grid. This is not: necessarily equal to the previous one, because for speed reasons (looking for points in expensive for N large) sometimes one point might be recomputed twice

adapted.N : the current number of dimensions considered for exploration adapted.private : a structure contained more detailed information on the status of the adaptive algorithm, that is needed

to resume the computation. In particular, the needed data structure consists of: private.G : the multi-idx set used to build the sparse grid. This is called

I in Gerstner-Griebel paper.

private.Ithe set of explored indices; G = I plus the neighbours of I.
This is called O in Gerstner-Griebel paper. I is sorted lexicographically

private.G_logsame as private.G, but sorted in the order with which indices
are added to G insted of lexicographic

private.I_logsame as private.I, but sorted in the order with which indices
are added to I insted of lexicographic

private.coeff_G : the coefficients of the indices in private.G in the combination technique private.N_log : for each iteration, the value of N_curr private.idx=idx : the idx whose neighbour is the next to be explored private.maxprof : the corresponding profit private.profits : the set of idx whose profit has been computed. They have been added to the grid

but their neighbour is yet to be explored.

private.idx_binthe corresponding set of profits. This is called A in
Gerstner-Griebel paper

private.var_with_ptsvector of variables in which we have actually put points.
length(var_with_pts) + controls.var_buffer_size = N_curr

private.nb_pts_log : for each iteration, the current nb_pts private.Hr : for non-nested points, the list of points visited by the

algorithm. Empty for nested-points

private.f_on_Hr : for non-nested points, the evaluations of f on Hr. Empty for nested-points

sgpykit.main.compute_sobol_indices_from_sparse_grid(S, Sr, nodal_values, domain, flags)[source]

Compute Sobol indices from a sparse grid approximation of a scalar-valued function.

This function computes the Sobol indices of a scalar-valued function f:R^N -> R in two steps: 1) Converts the sparse grid approximation of that function into its equivalent Polynomial Chaos Expansion (PCE). 2) Performs algebraic manipulations of the PCE coefficients.

Parameters:

Sdict: Sparse grid structure.
Srdict: Reduced sparse grid structure.
nodal_valuesnumpy.ndarray: Nodal values of the function on the sparse grid.
domaindict: Domain information.
flagsdict: Flags for the conversion to modal.

Returns:

Sob_inumpy.ndarray: Principal Sobol indices.
Tot_Sob_inumpy.ndarray: Total Sobol indices.
Meanfloat: Expected value of the function.
Varfloat: Variance of the function.

sgpykit.main.convert_to_modal(S, Sr, nodal_values, domain, flags)[source]

Convert a sparse grid interpolant to a modal (spectral) expansion.

This function recasts a sparse grid interpolant as a sum of orthogonal polynomials, i.e., computes the spectral expansion of the interpolant.

Parameters:

Sarray_like

Sparse grid structure.

Srarray_like

Reduced sparse grid structure.

nodal_valuesarray_like

Values of the target function evaluated on the reduced sparse grid. If the function is vector-valued (F: R^N -> R^V), then nodal_values is a matrix of size V x M, where M is the number of points in Sr.knots. If the function is scalar-valued, nodal_values is a row vector.

domainarray_like

Domain specification for the sparse grid. The format depends on the polynomial family: - For ‘legendre’, ‘chebyshev’: 2xN matrix [a1, a2, …; b1, b2, …]

defining the hyper-rectangle bounds.

For ‘hermite’: 2xN matrix [mu1, mu2, …; sigma1, sigma2, …] defining normal distribution parameters.
For ‘laguerre’: 1xN matrix [lambda1, lambda2, …] defining exponential distribution parameters.
For ‘generalized laguerre’: 2xN matrix [alpha1, alpha2, …; beta1, beta2, …] defining Gamma distribution parameters.
For ‘jacobi’: 4xN matrix [alpha1, alpha2, …; beta1, beta2, …; a1, a2, …; b1, b2, …] defining Beta distribution parameters on intervals [a_n, b_n].
For mixed families: cell array of length N, where each cell contains the parameters for the corresponding polynomial family.

flagsstr or list of str

Polynomial family or families to use for the expansion. Supported values are: ‘legendre’, ‘chebyshev’, ‘hermite’, ‘laguerre’, ‘generalized laguerre’, ‘jacobi’. For mixed families, provide a list of strings.

Returns:

modal_coeffsndarray: Matrix of modal coefficients. Each row corresponds to a multi-index and has V components.
Kndarray: Matrix of multi-indices, one per row.

Notes

The function processes each tensor grid in the sparse grid, computes the modal coefficients for each grid, and aggregates the coefficients by summing those with identical multi-indices.

sgpykit.main.create_sparse_grid(N, w, knots, lev2knots, idxset=None, S2=None)[source]

Generate a sparse grid (and corresponding quadrature weights) as a linear combination of full tensor grids, employing formula (2.9) of [Nobile-Tempone-Webster, SINUM 46/5, pages 2309-2345].

Parameters:

Nint: Number of dimensions.
wint: Integer non-negative value defining the multiindex set.
knotscallable or list of callables: Function(s) to generate knots in each direction. If a single function, it is used in every direction. The header of knots_function is [x,w]=knots_function(m).
lev2knotscallable or list of callables: Function(s) defining the relation between level and number of knots. If a single function, it is used in every direction. The header of m_function is m=m_function(i).
idxsetcallable, optional: Function defining the multiindex set. Default is sum(i-1).
S2StructArray, optional: Another sparse grid to recycle tensor grids from.

Returns:

SStructArray

Structure containing the information on the sparse grid (vector of tensor grids). Each tensor grid S[j] is a structure with the following fields:

knotsndarray
Vector containing the tensor grid knots.

weightsndarray
Vector containing the corresponding weights.

sizeint
Size of the tensor grid, S[j].size = prod(m).

knots_per_dimlist
Cell array (N components), each component is the set of 1D knots used to build the tensor grid.

mndarray
The input vector m, m == lev2knots(idx), m(i)==length(S[j].knots_per_dim[i]).

coeffint
How many times the tensor grid appears in the sparse grid (with sign).

idxndarray
The multiidx vector corresponding to the current grid.

Cndarray

Multi-index set used to generate the sparse grid. It is sorted lexicographically and contains all multi-indices, even those whose coefficient in the combination technique is 0.

Notes

The reduced multi-idx set D with only the indices with non-zero coefficient can be obtained as: ` nb_idx = len(S) D = np.vstack(S.idx).T # D.shape is (N, nb_idx) ` and will be lexicographic sorted as well.

sgpykit.main.create_sparse_grid_add_multiidx(new_idx, S_in, I_in, coeff_in, knots, lev2knots)[source]

Produce a grid obtained by adding a single multi-idx to a previously existing grid.

Parameters:

new_idxarray_like: An index that must be admissible with respect to the index set I_in. This condition won’t be checked by the function.
S_inarray_like: A sparse grid to which new_idx should be added.
I_inarray_like: The index set that was used to create S_in, either implicitly (by defining the rule in CREATE_SPARSE_GRID) or explicitly (by using CREATE_SPARSE_GRID_MULTIIDX). Note that this cannot be produced as the union of the indices S_in.idx, because S_in contains only tensors whose coefficient in the combination technique is non-zero, while here we need them all.
coeff_inarray_like: The vector of coefficients of the combination technique obtained on I_in. Note that this vector of coefficients MUST include also zeros if a row of I_in has coeff zero in the combination technique. Therefore, coeff_in cannot be taken as [S_in.coeff]. Instead, use coeff_in = COMBINATION_TECHNIQUE(I_in).
knotsarray_like or callable: The usual argument to specify knots.
lev2knotsarray_like or callable: The usual argument to specify lev2knots.

Returns:

Sarray_like: The new sparse grid.
Iarray_like: The new multiidx set, i.e., I = sortrows([I_in; new_idx]).
coeffarray_like: The updated vector of coefficients of the combination technique.

sgpykit.main.create_sparse_grid_multiidx_set(C, knots, lev2knots, S2=None, base=0)[source]

Create a sparse grid starting from a multiindex-set rather than from a rule IDXSET(I) <= W.

This function produces a sparse grid starting from a multiindex-set rather than from a rule IDXSET(I) <= W.

Parameters:

Carray_like: The multiindex set C. C must be in lexicographic order and admissible.
knotscallable or list of callable: Function(s) to generate knots for each dimension.
lev2knotscallable: Function to map levels to number of knots.
S2dict, optional: Another sparse grid. If provided, the function tries to recycle tensor grids from S2 to build those of S instead of recomputing them. This can be helpful whenever sequences of sparse grids are generated. Note that NO check will be performed whether S2 was generated with the same lev2knots as the one given as input. S2 can also be empty, S2=[].
baseint, optional: if matrix C is using 1-based indexing

Returns:

Sdict: The sparse grid structure.
Carray_like: The multiindex set C, possibly modified.

See also

check_set_admissibility: For admissibility check.
create_sparse_grid: For further information on KNOTS, LEV2KNOTS, and on the sparse grid data structure S.

sgpykit.main.create_sparse_grid_quick_preset(N, w)[source]

Create a sparse grid with a quick preset configuration.

This function creates a “vanilla” sparse grid using Clenshaw-Curtis points in the interval [-1, 1] with the lev2knots_doubling rule and a multi-index set defined by sum(ii) <= w. It is a shortcut for the following operations:

knots = lambda n: knots_CC(n, -1, 1)
S, C = create_sparse_grid(N, w, knots, lev2knots_doubling)
Sr = reduce_sparse_grid(S)

Parameters:

Nint: Number of dimensions.
wint: Level of the sparse grid.

Returns:

Sstruct: The full sparse grid structure.
Srstruct: The reduced sparse grid structure.
Cndarray: Coefficient matrix for the combination technique.

Notes

This function is equivalent to using the ‘SM’ (Smolyak) rule in define_functions_for_rule.

sgpykit.main.derive_sparse_grid(S, Sr, values_on_grid, domain, eval_points, h=None)[source]

Compute derivatives (gradients) of a scalar-valued function f: R^N -> R by centered finite differences formulas applied to the sparse grid approximation of f. The gradients can be computed at several points simultaneously.

Parameters:

Sstruct: Sparse grid struct.
Srstruct: Reduced version of S.
values_on_gridndarray: Values of the interpolated function on Sr (VALUES_ON_GRID is a vector, because the function is scalar-valued).
domainndarray: 2xN matrix = [a1, a2, a3, …; b1, b2, b3, …] defining the lower and upper bound of the hyper-rectangle on which the sparse grid is defined.
eval_pointsndarray: Points where the derivative must be evaluated. It is a matrix with points stored as columns, following the convention of the package.
hndarray or float, optional: Finite differences increment. H can be a scalar or a vector, in which case the n-th entry will be used as increment to approximate the n-th component of the gradient. If None, the increment size is chosen according to the length of each interval [an bn] as h_n = (b_n - a_n)/1E5.

Returns:

gradsndarray: Computed values of the derivatives (gradients). The gradients in each point are stored as columns of GRADS, i.e., size(GRADS) = N x size(EVAL_POINTS,2).

sgpykit.main.evaluate_on_sparse_grid(f, S, Sr, evals_old=None, S_old=None, Sr_old=None, tol=1e-14)[source]

Evaluate a function on a sparse grid, possibly recycling previous calls.

Several input combinations are possible, every input from the 3rd on can be set to [].

Parameters:

fcallable: Function to evaluate. Takes a column vector point and returns a scalar or column vector.
Ssparse grid or None: Sparse grid structure. If None, Sr must be provided.
Srreduced sparse grid or None: Reduced sparse grid structure. Must be provided if S is None.
evals_oldndarray, optional: Matrix storing previous evaluations of f on Sr_old.
S_oldsparse grid or None, optional: Previous sparse grid structure.
Sr_oldreduced sparse grid or matrix, optional: Previous reduced sparse grid or matrix of points.
tolfloat, optional: Tolerance for testing point equality (default 1e-14).

Returns:

f_evalndarray: Matrix storing evaluations of f on Sr.
new_pointsndarray: Points where f has been evaluated (new points w.r.t. previous grid).
tocomp_listndarray: Position of new_points in Sr.knots.
discard_pointsndarray: Points of Sr_old that have been discarded.
discard_listndarray: Index vector s.t. Sr_old.knots[:, discard_list] == discard_points.

Notes

Possible calls: - res, *_ = evaluate_on_sparse_grid(f, S=None, Sr) - res, *_ = evaluate_on_sparse_grid(f, S, Sr, evals_old, S_old, Sr_old) - res, *_ = evaluate_on_sparse_grid(f, S, Sr) - res, *_ = evaluate_on_sparse_grid(f, S, Sr, evals_old, None, Sr_old) - res, *_ = evaluate_on_sparse_grid(f, S, Sr, evals_old, S_old, Sr_old, tol)

API is like in evaluate_on_sparse_grid(). S or Sr have to be provided explicitly, opposed to SGMK.

sgpykit.main.export_sparse_grid_to_file(Sr, filename='points.dat', with_weights=False)[source]

Save knots of a reduced sparse grid to an ASCII file.

The first line of the file shows the number of points in the grid and their dimension. Then, points are stored as lines. If with_weights is True, the corresponding weight is added as the last entry of the row.

Parameters:

Srobject: Reduced sparse grid object containing knots and weights.
filenamestr, optional: Name of the file to save the points to. Default is ‘points.dat’.
with_weightsbool, optional: If True, weights are saved alongside the points. Default is False.

Notes

If with_weights is False, the file format is:
num_points dimension coord1(P1) coord2(P1) … coordN(P1) coord1(P2) coord2(P2) … coordN(P2) …
If with_weights is True, the file format is:
num_points dimension coord1(P1) coord2(P1) … coordN(P1) weight(P1) coord1(P2) coord2(P2) … coordN(P2) weight(P2) …

sgpykit.main.hessian_sparse_grid(S, Sr, values_on_grid, domain, eval_point, h=None)[source]

Compute the Hessian of a scalar-valued function f: R^N -> R by finite differences applied to the sparse grids approximation of f.

This function can only be evaluated at a single point, unlike derive_sparse_grid.

Parameters:

Sstruct: Sparse grid struct.
Srstruct: Reduced version of S.
values_on_gridndarray: Values of the interpolated function on Sr.
domainndarray: 2xN matrix = [a1, a2, a3, …; b1, b2, b3, …] defining the lower and upper bound of the hyper-rectangle on which the sparse grid is defined.
eval_pointndarray: Column vector point where the derivative must be evaluated.
hscalar or ndarray, optional: Finite differences increment size. If None, the increment size is chosen according to the length of each interval [an bn] as h_n = (b_n - a_n)/1E5.

Returns:

Hessianndarray: Hessian matrix, H(i,j) = partial_{y_i} partial_{y_j} sparse grid.

sgpykit.main.interpolate_on_sparse_grid(S, Sr, function_on_grid, non_grid_points)[source]

Interpolate a function on a sparse grid.

Evaluates the sparse grid polynomial approximation (surrogate model) on a generic point of the parameters space.

Parameters:

Slist of Struct or Struct: Sparse grid structure(s). If a single Struct is provided, it is converted to a list.
SrStruct: Reduced version of the sparse grid S.
function_on_gridndarray: Matrix containing the evaluation of the function on the points of Sr. Dimensions: V x number_of_points_in_the_sparse_grid.
non_grid_pointsndarray: Set of points where the sparse grid polynomial approximation is to be evaluated. Each column is a different point.

Returns:

ndarray: Matrix containing the evaluation of the vector-valued function in each of the non_grid_points. Dimensions: V x number_of_non_grid_point.

sgpykit.main.plot3_sparse_grid(S, dims, *args, **kwargs)[source]

Plot a sparse grid in 3D.

Parameters:

SStructArray: Sparse grid object to plot. S is a sparse grid in 3D. S can be either reduced or not. S can also be a tensor grid.
dimslist: List of dimensions to plot. It plots the components d1, d2, d3 of the points in S if S is more than 3D. If empty, defaults to [0, 1, 2].
*argstuple: Additional positional arguments for plot styling.
**kwargsdict: Additional keyword arguments for plot styling.

Returns:

hobject: Handle to the plot.

sgpykit.main.plot_sparse_grid(ax, S, dims, *args, **kwargs)[source]

Plot object based on the provided structure S and dimensions.

Plot a sparse grid in 2D or higher dimensions by selecting the specified dimensions. The sparse grid can be either reduced or not, or even a tensor grid.

Parameters:

axmatplotlib.axes.Axes: The axes on which to plot.
SStructArray: An object containing ‘knots’ as an attribute. Typically a sparse grid structure.
dimslist or None: Dimensions to plot (default is [1, 2]). If S is more than 2D, this specifies which dimensions to plot.
*argstuple: Additional positional arguments for plotting, such as format strings.
**kwargsdict: Additional keyword arguments for plotting, such as ‘color’, ‘marker’, etc.

Returns:

hmatplotlib.lines.Line2D: The handle to the plotted object.

sgpykit.main.plot_sparse_grids_interpolant(ax, S, Sr, domain, f_values, *args, **kwargs)[source]

Plot the sparse grid interpolant of a function.

Different plots are produced depending on the number of dimensions of the sparse grid: - If N==2, a surf plot will be generated. - If N==3, a number of contourfs (i.e., flat surfaces colored according to the value of the interpolant)

will be stacked over the same axes.

If N>3, a number of bidimensional cuts will be considered, and for each of them a surf will be generated. In other words, all variables but two will be frozen to their average value and the resulting two-dimensional plot will be produced.

Parameters:

axmatplotlib.axes.Axes: The axes on which to plot.
Sobject: The sparse grid.
Srobject: The reduced sparse grid.
domainnumpy.ndarray: A matrix 2 x N, describing the domain where the sparse grids should be plotted. The first row contains the lower bounds, the second row the upper bounds.
f_valuesnumpy.ndarray: The function values at the sparse grid points.
*argstuple: Additional arguments to control the behavior of the plots.
**kwargsdict: Additional keyword arguments to control the behavior of the plots.

Returns:

matplotlib.collections.PolyCollection or matplotlib.collections.LineCollection: The plot handle.

Notes

Additional inputs can be passed to control the behavior of the plots. Any combination of these optional inputs is allowed: - ‘with_f_values’: Adds dots with the values of the sparse grids interpolant to the plots above (case N=2 and N>3).

For N==3, adds the sparse grid points in the 3D plot.

‘nb_plot_pts’: Sets the number of points used in each direction for the surf/contourf plots (default 20).
‘nb_contourfs’: Sets the number of contourfs in the vertical direction for the case N=3 (default 5).
‘nb_contourf_lines’: Sets the number of contourf lines (default 10).
‘two_dim_cuts’: Specifies the couples of variables to consider for the two-dimensional cuts when N>3. C is a vector with 2*k components denoting the directions of the cuts. For instance, the default value is C = [0, 1, 2, …] and produces cut plots for (y1,y2) (y3,y4), (y5,y6).

sgpykit.main.quadrature_on_sparse_grid(f, S, Sr, evals_old=None, S_old=None, Sr_old=None, tol=None)[source]

Compute the integral of a function using a sparse grid.

This function behaves similarly to evaluate_on_sparse_grid, but returns the approximated integral of the function. See evaluate_on_sparse_grid for more information on inputs.

Parameters:

fcallable or array_like: Function handle or vector/matrix containing the evaluations of the function.
Sobject, optional: Sparse grid object.
Srobject: Reduced sparse grid object.
evals_oldarray_like, optional: Previous evaluations of the function.
S_oldobject, optional: Previous sparse grid object.
Sr_oldobject, optional: Previous reduced sparse grid object.
tolfloat, optional: Tolerance for incremental updates.

Returns:

resfloat: Approximated integral of the function.
evalsarray_like: Evaluations of the function over the points of the sparse grid.

Notes

Possible calls: - res, evals = quadrature_on_sparse_grid(f, S=None, Sr) - res, evals = quadrature_on_sparse_grid(f, S, Sr, evals_old, S_old, Sr_old) - res, evals = quadrature_on_sparse_grid(f, S, Sr) - res, evals = quadrature_on_sparse_grid(f, S, Sr, evals_old, None, Sr_old) - res, evals = quadrature_on_sparse_grid(f, S, Sr, evals_old, S_old, Sr_old, tol)

API is like in evaluate_on_sparse_grid(). S or Sr have to be provided explicitly, opposed to SGMK.

sgpykit.main.reduce_sparse_grid(S, tol=1e-14)[source]

Reduce a sparse grid by removing duplicate points within a given tolerance.

Given a sparse grid stored as a list of tensor grids, this function generates a reduced structure containing only non-repeated points and corresponding weights.

Parameters:

Sdict or object: Sparse grid structure with fields: - knots: list of tensor grid knots - weights: list of corresponding weights
tolfloat, optional: Tolerance to identify coincident points (default is 1e-14)

Returns:

Srdict or object: Reduced sparse grid structure with fields: - knots: list of non-repeated knots - weights: list of corresponding weights - size: number of non-repeated knots - n: forward map from [S.knots] to Sr.knots - m: inverse map from Sr.knots to [S.knots]