SGMethods package

In tis section, we document the SGMethods package. We give documentation for all (or most of) the modules, classes, and functions in the package.

sgmethods.compute_error_indicators module

Module with functions to compute the a-posteriori error indicators for adaptive enlargement of the multi-index set (equivalently, enlargement of the sparse grid or enrichment of the polynomial approximation space).

Coming soon: A-posteriori error estimator from Guignard, Nobile, SIAM JNA (2018).

sgmethods.compute_error_indicators.compute_GG_indicators_RM(old_RM, old_estimators, mid_set, knots, lev2knots, f, SG, u_on_SG, yy_rnd, L2_err_param, u_interp, TP_inteporlant=<class 'sgmethods.tp_inteprolants.TPPwLinearInterpolator'>, n_parallel=1)

Compute refinement error indicators for adaptive sparse grid interpolation as in

Gerstner, T., Griebel, M. Dimension-Adaptive Tensor-Product Quadrature. Computing 71, 65-87 (2003). https://doi.org/10.1007/s00607-003-0015-5

The error indicator corresponding to \(\nu\in\mathcal{R}\mathcal{M}_{\Lambda}\) is given by:

\[\left\vert \Delta_{\nu} u\right\vert_{L^2_{\mu}(\Gamma)},\]

where \(\Delta_{\nu} u\) is the hierarchical surplus operator, \(u:\Gamma\rightarrow V\) is the function we are interpolating, and \(\mu\) is the Gaussian measure. We compute the indicators only on the reduced margin of the current multi-index set. NB Recall that their sum (even on the whole margin) is NOT an error estimator. This function recycles the previously computed values if the margin did not change in a neighbourhood.

Parameters:

  • old_RM (numpy.ndarray[int]) – Old reduced margin. It is a 2D array where each row is a multi-index.

  • old_estimators (numpy.ndarray[float]) – Error indicators computed on the old reduced margin.

  • mid_set (MidSet) – Current multi-index set.

  • knots (Callable[[int], numpy.ndarray[float]]) – Returns the nodes vector of input length.

  • lev2knots (Callable[[int], int]) – Given a level >=0, returns a corresponding number number >0 of nodes.

  • f (Callable[[numpy.ndarray[float]], numpy.ndarray[float]) – Function to interpolate.

  • SG (numpy.ndarray[float]) – Current sparse grid. A 2D array.

  • u_on_SG (numpy.ndarray[float]) – Values F on current sparse grid. 2D array, each row is a value of F on a sparse grid node.

  • yy_rnd (numpy.ndarray[float]) – Random parameters in parametric domain. A N+1-dimensional array where each row is a parameter vector (dimension N).

  • L2_err_param (Callable[[numpy.ndarray[float], numpy.ndarray[float]], float]) – Compute the \(L^2(\Gamma)\) distance of the given functions using Monte Carlo quadrature. The functions are given through their value on the same random points in \(\Gamma\).

  • u_interp (numpy.ndarray[float]) – Values of the function to interpolate on the current sparse grid. 2D array, each row is the value in the corresponding row of SG.

  • TP_inteporlant (Appropriate tensor product interpolant class, optional) – Desired interpolation method as in tp_interpolants. Defaults to piecewise linear.

  • n_parallel (int, optional) – Number of parallel computations. Default 1.

Returns:

1D array. Each entry is the error indicator corresponding to a multi-index in the reduced margin.

Return type:

numpy.ndarray[float]

sgmethods.mid_set module

Module with a class for adaptive multi-index sets.

class sgmethods.mid_set.MidSet(max_n=inf, track_reduced_margin=False)

Bases: object

Multi-index set that can be grown by adding multi-indices from the margin. Can be used for adaptivity, where one does not know a-priori the next best multi-index set, or for a priori profit-based approximation when the definition of the profit is not monotonic.

The multi-index set is initialized to \(\{0\}\).

We can add multi-indices using either enlarge_from_margin(), which adds a multi-index from the margin and does not increasing the dimensionality, or increase_dim(), which adds the first coordinate unit vector outside of the support of the multi-index set. If the maximum dimension maxN is reached, then the dimension-adaptivity is ignored.

enlarge_from_margin(idx_margin)

Add the margin multi-index in position idx_margin to the multi-index set. A recursive function keeps the multi-index set downward-closed. This function does not increase the dimensionality of the multi-index (see increase_dim() for that).

Parameters:

idx_margin (int) – Index inside margin of the multi-index to be added to the multi-index set.

Returns:

None

get_cardinality_margin()

Returns the number of multi_indices in the margin of the multi-index set.

Parameters:

None –

Returns:

Number of multi-indices in the margin

Return type:

int

get_cardinality_mid_set()

Returns the number of multi-indices in the multi-index set.

Parameters:

None –

Returns:

Number of multi-indices in the multi-index set.

Return type:

int

get_cardinality_reduced_margin()

Returns the number of multi_indices in the reduced margin of the multi-index set.

Parameters:

None –

Returns:

Number of multi-indices in the reduced margin

Return type:

int

get_dim()

Returns the dimension (i.e. number of coordinates with a nonzero index) of the multi-index set.

Returns:

Dimension of the multi index set.

Return type:

int

increase_dim()

Increase the dimensionality of the multi-index set by adding the first coordinate unit vector outside of the support of the multi-index set. If the maximum dimension maxN is reached, then the dimension-adaptivity is not considered.

Parameters:

None –

Returns:

None.

sgmethods.multi_index_sets module

A collection of functions to generate downward-closed multi-index sets.

sgmethods.multi_index_sets.aniso_smolyak_mid_set(w, N, a)

Generate the anisotropic Smolyak multi-index set:

\[\Lambda = \{\nu \in \mathbb{N}_0^N : \sum_{i=1}^N a_i\nu_i \leq w \}.\]

where we denoted \(\Lambda\) the multi-index set, and \(w, N, a\) the input w, N, a. \(a=(a_i)_{i\in\mathbb{N}}\) is a sequence of positive real numbers that determines the anisotropy of the multi-index set.

Parameters:

  • w (float) – Sizing parameter of the multi-index set.

  • N (int) – Number of dimensions.

  • a (Callable[[int], numpy.array[float]]) – The resulting \(a(N)\) gives the anisotropy vector of length \(N\). The anisotropy vector is a vector of positive floats.

Returns:

Multi-index set (each row is a multi-index of length N).

Return type:

numpy.ndarray[int]

sgmethods.multi_index_sets.aniso_smolyak_mid_set_free_N(w, a)

Generate the anisotropic Smolyak multi-index set:

\[\Lambda = \{\nu\in\mathbb{N}_0^N : \sum_{i=1}^N a_i\nu_i \leq w\},\]

where we denoted \(\Lambda\) the multi-index set, and \(w, a\) the input w, a. \(a=(a_i)_{i\in\mathbb{N}}\) is a sequence of positive real numbers that determines the anisotropy of the multi-index set. The number of dimensions of the multi-index set is determined as:

\[N = \min \{ N\in\mathbb{N} : \forall \nu\in\Lambda,\ \sum_{n=1}^N a_n\nu_n \leq w \}.\]

Thus, the user must ensure that this quantity is finite!

Parameters:

  • w (float) – Sizing parameter of the multi-index set.

  • a (Callable[[int], numpy.ndarray[float]]) – a(N) gives the anisotropy vector of length \(N\). The anisotropy vector is a vector of positive floats.

Returns:

Multi-index set (each row is a multi-index).

Return type:

numpy.ndarray[int]

sgmethods.multi_index_sets.compute_mid_set_fast(profit_fun, minimum_profit, N)

Computes recursively over N dimensions the multi-index set such that the profit profit_fun is above the threshold minimum_profit for all its multi-indices. The multi-index set is defined as:

\[\Lambda = \{\nu\in\mathbb{N}_0^{N}:\mathcal{P}_{\nu}\geq P_{\textrm{min}}\}\]

where we denoted \(\Lambda\) the multi-index set, \(\mathcal{P}_{\nu}\) the profit function to apply to the multi-indices, and \(P_{\textrm{min}}>0\) the minimum profit threshold. The profit function should be monotone (decreasing) with respect to to the multi-indices to obtain a downward-closed multi-index set.

Parameters:

  • profit_fun (Callable[[numpy.ndarray[int]], float]) – Compute profit (positive float) for a multi-index.

  • minimum_profit (float) – minimum profit threshold for multi-index to be included (int): number of dimensions.

  • N (int) – Number of dimensions.

Returns:

Multi-index set (each row is a multi-index of length N).

Return type:

numpy.ndarray[int]

sgmethods.multi_index_sets.mid_set_profit_free_dim(profit_fun, minimum_profit)

Compute a multi-index set based on a profit function profit_fun and a minimum profit threshold profit_fun as:

\[\Lambda = \{\nu\in\mathbb{N}_0^{\infty}: \mathcal{P}_{\nu}\geq P_{\textrm{min}}\},\]

where we denoted \(\Lambda\) the multi-index set, \(\mathcal{P}_{\nu}\) the profit function applied to a multi-index, and \(P_{\textrm{min}}>0\) the minimum profit threshold. The profit should be monotone (decreasing) with respect to to the dimensions (otherwise it is not possible to determine it) and the multi-indices (to obtain a downward-closed multi-index set). This condition cannot be checked by the library, the user must ensure it. The present function also determines the maximum dimension possible, then calls the compute_mid_set_fast() with fixed dimension.

Parameters:

  • profit_fun (Callable[[numpy.ndarray[int]], float]) – Get profit (positive) for a multi-index. This function is required to vanish for increasing dimensions: \(\mathcal{P}_{\mathbf{e}_n} \rightarrow 0\) as \(n \rightarrow \infty\). Otherwise, the algorithm cannot terminate. This condition cannot be verified by the library and the user must ensure it.

  • minimum_profit (float) – minimum profit threshold for multi-index to be included

Returns:

Multi-index set (each row is a multi-index).

Return type:

numpy.ndarray[int]

sgmethods.multi_index_sets.smolyak_mid_set(w, N)

Generate the (isotropic) Smolyak multi-index set, i.e.

\[\Lambda = \{\nu \in \mathbb{N}_0^N : \sum_{i=1}^N \nu_i \leq w\},\]

where we denoted \(\Lambda\) the multi-index set, and \(w, N\) the input w, N.

Parameters:

  • w (float) – Sizing parameter of the multi-index set.

  • N (int) – Number of dimensions.

Returns:

Multi-index set (each row is a multi-index of length N).

Return type:

numpy.ndarray[int]

sgmethods.multi_index_sets.tensor_product_mid_set(w, N)

Generate a tensor product multi-index set, i.e.

\[\Lambda = \{\nu\in\mathbb{N}_0^N : \sum_{n=1}^N \nu_n \leq w\},\]

where we denoted \(\Lambda\) the multi-index set, and \(w, N\) the input w, N.

Parameters:

  • w (float) – Sizing parameter of the multi-index set.

  • N (int) – Number of dimensions.

Returns:

Multi-index set (each row is a multi-index of length N).

Return type:

numpy.ndarray[int]

sgmethods.multilevel_interpolant module

Module for a class that implements multilevel sparse grid interpolation. Multilevel methods combine a parametric and physical approximation method in a “smart” way in order to prodice a smaller error for the same number of degreesof freedom.

For a general introduction, see the classical paper on Mulilevel Monte Carlo:

Giles MB. Multilevel Monte Carlo methods. Acta Numerica. 2015;24:259-328. doi:10.1017/S096249291500001X, https://www.cambridge.org/core/journals/acta-numerica/article/abs/multilevel-monte-carlo-methods/C5AF9A57ED8FF8FDF08074C1071C5511

class sgmethods.multilevel_interpolant.MLInterpolant(interpolants_sequence)

Bases: object

Class for multilevel sparse grid interpolant as in

Teckentrup, Jantsch, Webster, Gunzburger. A Multilevel Stochastic Collocation Method for Partial Differential Equations with Random Input Data, SIAM JUQ, 2015 https://epubs.siam.org/doi/abs/10.1137/140969002.

We aim at approximating a function \(u:\Gamma\times D\rightarrow \mathbb{R}\), where \(\Gamma\subset\mathbb{R}^N\) is the parameter domain, \(D\subset\mathbb{R}^3\) is the physical domain (a time variable can also be included), and the codomain can be sobstituted by any Hilbert space.

The multilevel interpolant with \(K+1\) levels reads:

\[u_{\text{ML}}^K = \sum_{k=0}^{K} (I_{K-k}-I_{K-k-1}) u_k,\]

where \(I_{k}:\Gamma\rightarrow\mathbb{R}\) is a sparse grid interpolat of “resolution” \(k=0,\dots, K\), and \(u_k:\Gamma\times D\rightarrow\mathbb{R}\) is a space approximation (e.g. finite elements) with “resolution” \(k=0,\dots, K\).

get_ml_terms(yy, ml_samples_f)

Get the terms of the multi-level expansion split based on FE approximation.

Parameters:

  • yy (numpy.array[float]) – Each row is a parametric point where to evaluate each of the multilevel terms.

  • ml_samples_f (list[numpy.ndarray[float]]) – The k-th list entry contais data on \((I_{K-k}-I_{K-k-1})[u_k]\). Each row contains the finite element coordinates ina point in the level K-k sparse grid.

Returns:

List of 2D array. The k-th list entry contains the result for \((I_{K-k}-I_{K-k-1})[u_k]\). Each row contains the finite element coordinates ina point of yy.

Return type:

[list[numpy.ndarray[float]]]

interpolate(yy, ml_samples_f)

Use the interpolation operator to approximate the function f of which ML samples are given.

Parameters:
Returns:

2D array. Each row is the approximation of f in the parameter given by the corresponding row of yy.

Return type:

numpy.ndarray[double]

sample(f_approx)

Compute multilevel samples of the function f on appopriate sparse grids noeed to compute the multilevel interpolant.

Parameters:

f_approx (Callable[[numpy.ndarray[float], int], float]) – Function with 2 inputs: Parameters (2D numpy.ndarray[float], where each row is a parameter vector), and the finite element accuracy indexed by \(k=0,\dots, K\).

Returns:

List of the same length as self.interp_seq used to interpolate with interpolate(). The k-th list element is a 2D numpy.ndarray[float] containing the values of f_approx[k] (the level k finite element approximation) on the sparse grid of resolution K-k-1. Within this 2D array, each row corresponds to a point in the sparse grid at level K-k and consists of the coorodinates of the finite element approximation in that parameter.

Return type:

list[numpy.ndarray[float]]

total_cost(cost_kk)

Compute the cost of computing the multilevel interpolant.

Parameters:

cost_kk (numpy.ndarray[float]) – The k-th entry is the cost of computing 1 finite element sample at level k.

Returns:

Total cost based on number of sparse grid nodes, levels, and finite element computations.

Return type:

float

sgmethods.nodes_1d module

This module provides various functions to generate 1D nodes for interpolation and numerical integration. Dofferent nodes are usually suited to different collocation approximations.

sgmethods.nodes_1d.cc_nodes(n)

Generate n Clenshaw-Curtis nodes, i.e. extrema of (n-1)-th Chebyshev polynomial on [-1,1]. They read

\[x_i = -\cos\left(\frac{i\pi}{n-1}\right), \quad i=0,1,...,n-1.\]
Parameters:

n (int) – number of nodes >=1.

Returns:

n nodes; if n=1, return 0.

Return type:

numpy.ndarray[float]

sgmethods.nodes_1d.equispaced_interior_nodes(n)

Generates n equispaced nodes on (-1,1) with the first and last node in the interior at same distance from -1 and 1. They coincide with the n+2 equispaced nodes on [-1,1] with first and last node left out.

Parameters:

n (int) – number of nodes.

Returns:

n nodes.

Return type:

numpy.ndarray[float]

sgmethods.nodes_1d.equispaced_nodes(n)

Gnerate n equispaced nodes on [-1,1] with first and last node on boundary.

Parameters:

n (int) – number of nodes.

Returns:

n nodes.

Return type:

numpy.ndarray[float]

sgmethods.nodes_1d.hermite_nodes(n)

Hermite interpolation nodes (roots of the n-th Hermite polynomial).

Parameters:

n (int) – number of nodes >=1.

Returns:

n nodes; if n=1, return 0.

Return type:

numpy.ndarray[float]

sgmethods.nodes_1d.optimal_gaussian_nodes(n, p=2)

Generate n interpolation nodes that eventually cover \(\mathbb{R}\) and give optimal degree p+1 piecewise poynomial interpolation in the \(L^2_{\mu}(\mathbb{R})\) norm, where \(\mu\) denote the Gaussian density.

Parameters:

  • n (int) – number of nodes.

  • p (int, optional) – Piecewise polynomial interpolant degree + 1 (i.e. the number of nodes needed to determine a piecewise polynomial of degree p in each interval). Defaults to 2. n must be odd because nodes are symmetric around 0 and one node is 0.

Returns:

n interpolation nodes.

Return type:

numpy.ndarray[float]

sgmethods.nodes_1d.unbounded_nodes_nested(n, p=2)

Use optimal_gaussian_nodes() to generated n nested nodes that are, in the limit, dense in \(\mathbb{R}\) and give optimal degree p+1 piecewise poynomial interpolation in the \(L^2_{\mu}(\mathbb{R})\) norm, where \(\mu\) denote the Gaussian density. n must be odd because one node is 0 and nodes are symmetric around 0. n must have the right value for the sequence of nodes to be nested, i.e. \(n = 2^{i+1}-1\) for some \(i\in\mathbb{N}_0\).

Parameters:

  • n (int) – number of nodes.

  • p (int, optional) – Piecewise polynomial interpolant degree + 1 (i.e. the number of nodes needed to determine a piecewise polynomial of degree p in each interval). Defaults to 2.

Returns:

n interpolation nodes.

Return type:

numpy.ndarray[float]

sgmethods.nodes_tp module

This module provides the function to generate tensor product nodes from 1D nodes.

sgmethods.nodes_tp.tp_knots(scalar_knots, num_knots_dir)

Generate tensor product nodes given 1D nodes family and number of nodes in each direction.

Parameters:
Returns:

Tensor product nodes. Each tuple entry contains the nodes in one direction.

Return type:

tuple[numpy.ndarray[float]]

sgmethods.parametric_expansions_brownian_motion module

This module provides functions for parametric expansions of the Wiener process.

Functions:
param_LC_Brownian_motion(tt, yy, T):

Constructs the Wiener process using the classical Levy-Ciesielski construction.

param_KL_Brownian_motion(tt, yy):

Constructs the Wiener process using the Karhunen-Loeve expansion.

sgmethods.parametric_expansions_brownian_motion.param_KL_Brownian_motion(tt, yy)

The Karhunen_Loeve expansion of th Brownian motion. Can be computed exactly.

Parameters:

  • tt (numpy.ndarray[float]) – Discrete times of evaluation in [0,1].

  • yy (numpy.ndarray[float]) – 2D array. Each row is a parameter vector of the expansion. Each component is a real numbers that replace i.i.d. standard Gaussians.

Returns:

Each row gives the samples of the Wiener process on tt for the corresponding row (parameter vector) in yy.

Return type:

numpy.ndarray[float]

sgmethods.parametric_expansions_brownian_motion.param_LC_Brownian_motion(tt, yy, T)

The classical Lèvy-Ciesielsky construction of the Wiener process is used as a parametric expansion as follows:

\[W(y, t) = y_0 \eta_{0,1}(t) + \sum_{l=1}^L \sum_{j=1}^{2^{l-1}} y_{l,j} \eta_{l,j}(t),\]

where \(\eta_{l,j}\) is the Faber-Schauder basis on \([0,T]\) (i.e. a wavelet basis of hat functions) and \(y = (y_{l,j})_{l,j}\) is a sequence of real numbers that replace i.i.d. standard Gaussians.

Parameters:

  • tt (numpy.ndarray[float]) – Discret times of evaluation in \([0,T]\).

  • yy (numpy.ndarray[float]) – The parameter vectors of the expansion. Each row consists of the scalar components (each in \(\mathbb{R}\)) of a parametric vector.

  • T (float) – Final (positive) time of approximation.

Returns:

Each row gives the approximation of the function in one parametric poin in yy through its values at the discrete sample times in tt.

Return type:

numpy.ndarray[float]

sgmethods.sparse_grid_interpolant module

The core module of SGmethods. It contains the implementation of sparse grid interpolation in a class.

class sgmethods.sparse_grid_interpolant.SGInterpolant(mid_set, knots, lev2knots, tp_interpolant=<class 'sgmethods.tp_inteprolants.TPPwLinearInterpolator'>, n_parallel=1, verbose=True)

Bases: object

Sparse grid interpolant class. It stores all relevant information to define it like multi-index set, 1D notes etc. It automatically computes the sparse grid and inclusion-exclusion coefficients upon initialization. It allows to interpolate high dimensinal functions on the sparse grid.

interpolate(x_new, f_on_SG)

Evaluate the interpolant on new parametric points.

Parameters:

  • x_new (numpy.ndarray[float]) – 2D array. New parametric points where to interpolate. Each row is a point.

  • f_on_SG (numpy.ndarray[float]) – 2D array. Values of function on the sparse grid. Each row is a value corresponding to a parametric point in the sparse grid self.SG.

Returns:

2D array. Values of the interpolant of f on x_new. Each row is a value corresponds to the parameter stored in the same row of x_new.

Return type:

numpy.ndarray[float]

sample_on_SG(f, dim_f=None, old_xx=None, old_samples=None)

Sample a given function on the sparse grid.

Optionally recycle previous samples stored e.g. from a previous interpolation. The class also takes care automatically of the case in which the sparse grid has increased dimension (number of approcimated scaar parameters). The class first checks whether or not there is anything to recycle, then it sample new values.

NB The class assumes f takes as input a numpy.ndarray[float] of parameters. The 1st output is assumed to be a numpy.ndarray[float] or a float (in this case, it is transformed into a 1-entry numpy.ndarray[float]).

Parameters:

  • f (Callable[[numpy.ndarray[float]], numpy.ndarray[float]]) – The function to interpolate.

  • dim_f (int, optional) – Dimensinoality codomain of f. Defaults to None.

  • old_xx (numpy.ndarray[float], optional) – 2D array. Each row is a parametric point Defaults to None.

  • old_samples (numpy.ndarray[float], optional) – 2D array. Each row corresponds to a row of old_xx. Defaults to None.

Returns:

2D array. Each row is the value of f on a sparse grid (self.SG) point.

Return type:

numpy.ndarray[float]

setup_SG()

Computes and saves in class attributes some important quantities: SG (sparse grid), num_nodes, map_tp_to_SG based on active_tp_nodes_list.

setup_interpolant()

Computes and saves in class attributes some important quantities: combination coefficients, active multi-indicies, active TP dimensions, active TP nodes, map TP to SG nodes based on the multi-index set, the nodes and the level-to-knot function.

Parameters:

None –

Returns:

None

sgmethods.tp_inteprolants module

Examples of tenosr products interpolatio operators. They will be wrapped in the interpolant wrapper class and used for sparse grid interpolation through the inclusion-exclusion formula. The user can implement new interpolation operators ina analogous classes. IMPORTANT: The input and output must be compatible with the wrapper, so same as the examples in this file! In particular the following methods are needed: - __init__(self, nodes_tuple, f_on_nodes); - __call__(self, x_new): Input is numpy.ndarray[float] (each row is a parameter vector on which to interpolate), output is numpy.ndarray[float] (each row is the value of the interpolant on the corresponding parameter in x_new).

class sgmethods.tp_inteprolants.TPLagrangeInterpolator(nodes_tuple, f)

Bases: object

Given function samples on tensor product nodes, interpolates with tensor product Lagrange interpolation computed with barycentric interpolation formula, see Trefetten-Approximation Theory and Approximation Practice (2019)

class sgmethods.tp_inteprolants.TPPwCubicInterpolator(nodes_tuple, f)

Bases: object

Given function samples on tensor product nodes, interpolates with tensor product piecewise cubic interpolation.

class sgmethods.tp_inteprolants.TPPwLinearInterpolator(nodes_tuple, f_on_nodes)

Bases: object

Given function samples on tensor product nodes, interpolates with tensor product piecewise linear interpolation.

class sgmethods.tp_inteprolants.TPPwQuadraticInterpolator(nodes_tuple, f)

Bases: object

Given function samples on tensor product nodes, interpolates with tensor product piecewise quadratic interpolation.

sgmethods.tp_interpolant_wrapper module

Module with the implementation of the wrapper class of tensor product interpolants. Important to enforce the compatibility of user-defined interpolants.

class sgmethods.tp_interpolant_wrapper.TPInterpolatorWrapper(active_nodes_tuple, active_dims, f_on_nodes, tp_interpolant)

Bases: object

Wrapper for several tensor product interpolants.

If a direction has only 1 collocation node, it should be 0 and the 1-node approximation is constant in that direction. This is the case of any of the directions in which the multi-index does not have its support (no non-zero index).

The user can define new tensor prodcut interplants following the instructions written in the module tp_interpolants and look at the examples given there.

sgmethods.utils module

A collection of functions that carry out efficiently various basic tasks needed in the core functions or that may come in handy when using SGMethods in numerical aexperiments and applications.

sgmethods.utils.checkIfElement(mid, mid_set)

DEPRECATED. To be removed soon and substitute with find_mid().

sgmethods.utils.compute_level_lc(i)

Converts from linear to hierarchical indices for a wavelet expansion like the Levy-Ciesielski expansion of the Brownian motion. The change of indices is, for a linear index \(i\in\mathbb{N}_0\):

\[\begin{split}l &= \left\lfloor \log_2(i+1) \right\rfloor + 1 \\ j &= i - 2^{l-1}\end{split}\]
Parameters:

i (int) – Linear index. The fist basis function corresponds to i=0.

Returns:

Pair of hierarchical indices (level \(l\), index in level \(j\)). The levels start from 0, the index in level from 1.

Return type:

tuple[int, int]

sgmethods.utils.coord_unit_vector(dim, entry)

Returns a coordinate unit vector, i.e. one with all zeros except for one entry equal to one.

Parameters:

  • dim (int) – Length of the vector.

  • entry (int) – The unique entry equal to one.

Returns:

The desired coordinate unit vector.

Return type:

numpy.ndarray[int]

sgmethods.utils.find_idx_in_margin(margin, mid)

DEPRECATED. To be removed soon and substitute with find_mid().

sgmethods.utils.find_mid(mid, mid_set)

Find out whether a given multi-index appears in a multi index set and determine the position of its first occurrence.

Parameters:

  • mid (numpy.ndarray[int]) – The multi-index to search for.

  • mid_set (numpy.ndarray[int]) – The multi-index set. Must be a 2D array, each row is a multi-index.

Returns:

First value tells whether mid was found inside mid_set. Second value is the position (-1 if not found).

Return type:

tuple(bool,int)

sgmethods.utils.is_downward(mid_set)

Check if the multi-index set is downward-closed. This means that:

\[\forall\nu\in\Lambda, \forall n\in\mathbb{N},\nu_n=0\text{ or } \nu-e_n\in\Lambda,\]

where \(e_n\) is the n-th coordinate unit vector.

Returns:

True if the multi-index set is downward-closed, False otherwise.

Return type:

bool

sgmethods.utils.lexic_sort(mid_set)

Sort given multi-index set in lexicographic order. This means that one vector comes before another if the first non-equal entry (from the left) is smaller.

Parameters:

mid_set (numpy.ndarray[int]) – The multi-index set to test. Must be a 2D array, each row is a multi-index.

Returns:

The sorted multi-index set.

Return type:

bool

sgmethods.utils.ls_fit_loglog(x, y)

Fit samples of a function \(f:\mathbb{R}\rightarrow\mathbb{R}\) to a linear function in the least-squares sense in the log-log scale.

Parameters:
Returns:

The slope and offset of the linear fit in the log- log scale.

Return type:

tuple[float, float]

sgmethods.utils.mid_is_in_reduced_margin(mid, mid_set)

Check if mid belongs to the reduced margin of mid_set. This is the case if and only if

\[\forall n = 1,...,N, \ \nu_n = 0 \textrm{ or } \nu - e_n \in \Lambda,\]

where we denoted mid_set and mid by \(\Lambda\) and \(\nu\) respectively. Here \(e_n\) is the n-th coordinate unit vector, \(N\) is the dimension of mid_set and mid.

Parameters:
Returns:

True if mid is in the reduced margin of mid_set, False otherwise.

Return type:

bool

sgmethods.utils.plot_2d_midset(midset)

Plot the 2D multi-index set in the plane.

sgmethods.utils.rate(err, n_dofs)

Compute the rate of logarithmic convergence of a sequence.

Parameters:
Returns:

1D array containing rates of convergence computed between subsequence pars of n_dofs and err as

\[-\frac{\log(err_{i+1}/err_i)} {\log(n_{\text{dofs},i+1}/n_{\text{dofs},i})}.\]

Return type:

numpy.ndarray[float]