Simplex Stochastic Collocation tutorial¶
Standard Stochastic Collocation (SC) uses global polynomials to construct a surrogate model $\tilde{f}$ of some computational model $f$:
\begin{align} f(\boldsymbol\xi) = f(\xi_1,\cdots, \xi_d) \approx \tilde{f}(\xi_1,\cdots, \xi_d) := \sum_{j_1=1}^{m_{l_1}}\cdots\sum_{j_d = 1}^{m_{l_d}}f\left(\xi_{j_1}, \cdots, \xi_{j_d}\right) a_{j_1}(\xi_1)\otimes\cdots\otimes a_{j_d}(\xi_d) \end{align}Here $a_{j_i}(\xi_i)$ are (typically) 1D Lagrange interpolation polynomials, which span the entire domain of $\xi_i$ (i.e. they are global interpolation polynomials).
Challenges exist if the number of parameters $d$ is high, which is the well-known curse of dimensionality. The SC method can be made dimension-adaptive however, to achieve must better scalability. Click here for a practical EasyVVUQ tutorial on how to use the dimension-adaptive SC sampler. The maths behind the dimension-adaptivity are described here.
However, this tutorial is on the 2nd main challenge of SC (and PCE) methods, namely the assumed regularity of the code output. If discontinuities or high gradient are present in the stochastic input space $\boldsymbol\xi$, the quality of the surrogate $\tilde f$ deteriorates. To see this, run the example below which creates a SC surrogate for a 2D discontinuous function. We will assume you are familiar with the basic EasyVVUQ building blocks. If not, you can look at the basic tutorial.
Install EasyVVUQ¶
If you have not installed EasyVVUQ, uncomment the line below
#%pip install EasyVVUQ
Standard SC on a discontinuous function¶
import os
import chaospy as cp
import numpy as np
import easyvvuq as uq
import matplotlib.pyplot as plt
from easyvvuq.actions import CreateRunDirectory, Encode, Decode, ExecuteLocal, Actions
from matplotlib import cm
This is the function we consider, it is located in discont_model.py as well.
def f(x, n_xi = 2):
"""
Discontinuous test function
Parameters
----------
x : array
input parameters
Returns
-------
float
f(x)
"""
if x[1] <= -0.6 * x[0] + 0.8:
return np.sum(x[0:n_xi]) - 1
else:
return x[0] ** 3 + np.cos(np.sum(x[1:n_xi] ** 2)) + 1
This is what it looks like:
xx = np.linspace(0,1,20)
X1, X2 = np.meshgrid(xx, xx)
F = np.zeros([20,20])
for i in range(20):
for j in range(20):
F[i, j] = f(np.array([xx[j], xx[i]]))
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d', xlabel=r'$\xi_1$', ylabel=r'$\xi_2$')
surf = ax.plot_surface(X1, X2, F, linewidth=0, antialiased=False, cmap=cm.coolwarm)
plt.tight_layout()
This is a standard EasyVVUQ SC campaign, where we prescribe a standard uniform distribution for inputs $\xi_1$ and $\xi_2$:
def get_campaign(n_xi = 2):
##################################
# Define (total) parameter space #
##################################
# Define parameter space
params = {}
for i in range(5):
params['xi_%d' % (i + 1)] = {'type': 'float', 'default': 0.5}
params['n_xi'] = {'type': 'integer', 'default': 2}
###########################
# Set up a fresh campaign #
###########################
encoder = uq.encoders.GenericEncoder(
template_fname= os.path.abspath('discont_model.template'),
delimiter='$',
target_filename='input.csv')
execute = ExecuteLocal('{}/discont_model.py'.format(os.getcwd()))
output_columns = ["f"]
decoder = uq.decoders.SimpleCSV(
target_filename='output.csv',
output_columns=output_columns)
actions = Actions(CreateRunDirectory(root='/tmp', flatten=True), Encode(encoder), execute,
Decode(decoder))
campaign = uq.Campaign(name='standard_SC', work_dir='/tmp', params=params, actions=actions)
return campaign
# get the main Campaign object
campaign = get_campaign()
#######################
# Specify input space #
#######################
vary = {
"xi_1": cp.Uniform(0.0, 1.0),
"xi_2": cp.Uniform(0.0, 1.0)
}
##################
# Select sampler #
##################
sampler = uq.sampling.SCSampler(vary=vary, polynomial_order=6)
# associate the sampler with the campaign
campaign.set_sampler(sampler)
###############################
# execute the defined actions #
###############################
campaign.execute().collate()
# get EasyVVUQ data frame
data_frame = campaign.get_collation_result()
############################
# Post-processing analysis #
############################
analysis = uq.analysis.SCAnalysis(sampler=sampler, qoi_cols=["f"])
results = analysis.analyse(data_frame=data_frame)
results.describe()
| f | |
|---|---|
| 0 | |
| mean | 0.972725 |
| std | 1.294315 |
| var | 1.675251 |
The displayed output stats do not immediately hint at a problem. However, when we plot the surrogate $\tilde f$ we find:
##############################
# Plot the surrogate samples #
##############################
def plot_surrogate(analysis):
if len(vary) == 2:
#generate n_mc samples from the input distributions
n_mc = 20
xx = np.linspace(0, 1, n_mc)
X, Y = np.meshgrid(xx, xx)
F = np.zeros([n_mc, n_mc])
# evaluate the surrogate at these values
for i in range(n_mc):
for j in range(n_mc):
xi_mc = np.array([X[i, j], Y[i, j]])
F[i, j] = analysis.surrogate('f', xi_mc)
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d', xlabel=r'$\xi_1$', ylabel=r'$\xi_2$')
surf = ax.plot_surface(X, Y, F, linewidth=0, antialiased=False, cmap=cm.coolwarm)
plt.tight_layout()
elif len(vary) == 1:
n_mc = 100
xx = np.linspace(0, 1, n_mc)
fig = plt.figure()
ax = fig.add_subplot(111, xlabel=r'$\xi$', ylabel=r'$\tilde{f}$')
F = np.zeros([n_mc])
for i in range(n_mc):
F[i] = analysis.surrogate('f', xx[i])
ax.plot(xx, F, label=r'$\tilde{f}$')
ax.plot(analysis.sampler.tri.points, analysis.samples['f'], 'o',label=r'$f$')
plt.legend()
plt.tight_layout()
plot_surrogate(analysis)
Clearly, $\tilde f$ contains many "unphysical" oscillations. This is Runge phenomenon, caused by trying to model a function $f$ with highly localized & abrupt changes, by a function $\tilde f$ based on global & smooth polynomials. In other words, for this function our choice of interpolation polynomials (the Lagrange polynomials $a_{j_i}(x_i)$ shown above) is not good. A localized polynomial basis is required, that is able to adapt the local behavior of $f$.
Simplex stochastic collocation¶
The Simplex Stochastic Collocation (SSC) method differs from traditional collocation methods (SC/PCE), as it employs the Delaunay triangulation to discretize the probability space into simplex elements, rather than relying on the more common tensor product of one-dimensional abscissas. Using a multi-element technique has the advantage that local mesh adaptation can be performed, such that only regions of interest are refined. Since the SSC method discretizes the parameter space $\Xi$ into $n_e$ disjoint simplices $\Xi = \Xi_1 \cup \cdots \cup \Xi_{n_e}$, the local surrogate function $w_j$ (for inputs $\boldsymbol\xi\in\Xi_j$) is given by the expansion
\begin{align} w_j(\boldsymbol\xi) = \sum^{N_j}_{l = 0} c_{j,l}\Psi_{j, l}(\boldsymbol\xi). \end{align}The $N_j+1$ number of points needed for $n_\xi$-dimensional interpolation ($n_\xi$ is the number of random inputs $\xi_i$, $i=1,\cdots,n_\xi$) of polynomial order $p_j$ is given by
\begin{align} \label{eq:Np1_j} N_j + 1 = \frac{(n_\xi+p_j)!}{n_\xi!p_j!}. \end{align}The interpolation basis functions $\Psi_{j, l}$ are simple monomials, and the coefficients $c_{j, l}$ are found by solving a Vandermonde system, see the references below.
The SSC method automatically reduces the polynomial order if a stencil $S_j$ crosses a discontinuity. It achieves this by enforcing the so-called Local Extremum Conserving (LEC) limiter to all simplices $\Xi_j$ in all $S_j$. The LEC condition is given by
\begin{align} \boxed{ \min_{\boldsymbol\xi\in\Xi_j} w_j(\boldsymbol\xi) = \min{\bf v}_j \wedge \max_{\boldsymbol\xi\in\Xi_j} w_j(\boldsymbol\xi) = \max{\bf v}_j}, \end{align}where ${\bf v}_j = \{v_{k_{j,0}},\cdots,{v_{k_{j,n_\xi}}}\}$ are the code samples at the vertices of $\Xi_j$. If $w_j$ violates the LEC condition above in one of its $\Xi_j\in S_j$, the polynomial order $p_j$ of that stencil is reduced, usually by 1. Since polynomial overshoots occur when trying to interpolate a discontinuity, $p_j$ is reduced the most in discontinuous regions. Interpolating a function on a simplex with $p_j = 1$ and ${\bf v}_j$ located at its vertices always satisfies the LEC condition. This ensures that $w(\boldsymbol\xi)$ is able to represent discontinuities without suffering from the Runge phenomenon. In practice, the LEC condition is enforced for all $\Xi_j $ in all $S_j$ via random sampling of the $w_j$. If for a given $w_j$ the LEC condition is violated for any of the randomly placed samples $\boldsymbol\xi_j$, the polynomial order of the corresponding stencil is reduced.
Which $N_j+1$ points are used for $w_j$ is determined by the interpolation stencil $S_j$. The stencil can be constructed based on the nearest-neighbour principle. In this case the first $n_\xi+1$ points are the vertices $\{\boldsymbol\xi_{k_{j,0}},\cdots,\boldsymbol\xi_{k_{j,n_\xi}}\}$ of the simplex $\Xi_j$, which would suffice for $p_j = 1$. For higher degree interpolation, neighbouring points $\boldsymbol\xi_k$ are added based on their proximity to the center of simplex $\Xi_j$, i.e. based on their ranking according to
\begin{align} \label{eq:center_j} \Vert \boldsymbol\xi_k - \boldsymbol\xi_{\mathrm{center},j} \rVert_2, \end{align}where those $\boldsymbol\xi_k$ of the current simplex $\Xi_j$ are excluded. The simplex center $\boldsymbol\xi_{\mathrm{center},j}$ is defined as
\begin{align} \boldsymbol\xi_{\mathrm{center},j} := \frac{1}{n_\xi+1}\sum^{d}_{l = 0}\boldsymbol\xi_{k_{j,l}}. \end{align}The nearest neighbor stencil leads to a $p_j$ distribution that can increase quite slowly when moving away from a discontinuity. An example of this behavior can be found in the figure below, which shows the Delaunay triangulation with a discontinuity running through the domain (dotted line). An alternative to nearest-neighbor stencils are stencils created according to the Essentially Non-Oscillatory (ENO) principle [1]. The idea behind ENO stencils is to have higher degree interpolation stencils up to a thin layer of simplices containing the discontinuity. For a given simplex $\Xi_j$, its ENO stencil is created by locating all the nearest-neighbor stencils that contain $\Xi_j$, and subsequently selecting the one with the highest $p_j$. This leads to a Delaunay triangulation which captures the discontinuity better than its nearest-neighbor counterpart.
Like dimension-adaptive sampling, the SSC method iteratively refines the sampling plan. Which simplex element is refined is based on a geometrical refinement measure, given by;
\begin{align} \label{eq:ref_meas} \bar{e}_j := \bar{\Omega}_j\bar{\Xi}_j^{2O_j}. \end{align}It contains the probability and the volume of simplex $\Xi_j$, i.e.
\begin{align} \label{eq:prob_and_vol} \bar{\Omega}_j = \int\limits_{\Xi_j}f_{\boldsymbol\xi}(\boldsymbol\xi)d\boldsymbol\xi\;\;\;\mathrm{and}\;\;\; \bar{\Xi}_j = \int\limits_{\Xi_j}d\boldsymbol\xi, \end{align}where $\bar{\Xi} = \sum^{n_e}_{j=1}\bar{\Xi}_j$. The probability $\bar{\Omega}_j$ can be computed by Monte-Carlo sampling and $\bar{\Xi}_j$ via the relation
\begin{align} \label{eq:vol_j} \bar{\Xi}_j = \frac{1}{n_\xi!}\lvert\mathrm{det}\left(D\right)\rvert,\;\;\;D = \left[\boldsymbol\xi_{k_{j,1}}-\boldsymbol\xi_{k_{j,0}}\;\;\;\;\boldsymbol\xi_{k_{j,2}}-\boldsymbol\xi_{k_{j,0}} \;\;\;\;\cdots\;\;\;\; \boldsymbol\xi_{k_{j,n_\xi+1}}-\boldsymbol\xi_{k_{j,0}}\right] \in\mathbb{R}^{n_\xi \times n_\xi}. \end{align}Finally, the order of convergence $O_j$ is given by
\begin{align} \label{eq:order_conv} O_j = \frac{p_j+1}{n_\xi}. \end{align}The simplex with the highest $\bar{e}_j$ is selected for refinement. To ensure a good spread of the sample points, only randomly-selected points inside a simplex sub-element $\Xi_{sub_j}$ are allowed. The vertices of this sub-element are defined as
\begin{align} \label{eq:sub_simplex} \boldsymbol\xi_{sub_{j,l}} := \frac{1}{n_\xi}\sum^{n_\xi}_{\substack{l^*=0 \\ l^*\neq l}} \boldsymbol\xi_{k_{j,l^*}}, \end{align}see the figure below for a two-dimensional example. The SSC algorithm can be partially parallelized by selecting the $N<n_e$ simplices with the $N$ largest $\bar{e}_j$ for refinement.
Note that $\bar{e}_j$ is probabilistically weighted through $\Omega_j$ and that it assigns high value to those simplices with low $p_j$ since in general $\bar{\Xi}_j\ll1$. Effectively this means that $\bar{e}_j$ is a solution-dependent refinement measure which refines simplices near discontinuities since the SSC method automatically reduces the polynomial order (via LEC) if a stencil $S_j$ crosses a discontinuity.
This is the overall overview of the SSC method. For more information consult the articles below and the references therein.
References
[1] Witteveen, J. A. S., & Iaccarino, G. (2013). Simplex stochastic collocation with ENO-type stencil selection for robust uncertainty quantification. Journal of Computational Physics, 239, 1-21.
[2] Edeling, W. N., Dwight, R. P., & Cinnella, P. (2016). Simplex-stochastic collocation method with improved scalability. Journal of Computational Physics, 310, 301-328.
Current limitations¶
- We do not recommend running the SSC sampler using more than (approximately) 5 inputs variables.
- At this point, only scalar output variables are supported.
Both limitations will be addressed in a later release of EasyVVUQ.
Simplex Stochastic Collocation on a discontinuous function¶
We generate the EasyVVUQ campaign in the exact same way as in the SC case;
campaign = get_campaign()
Selecting the SSC sampler is also very similar to the SC case;
##################
# Select sampler #
##################
sampler = uq.sampling.SSCSampler(vary=vary, max_polynomial_order=4)
# Associate the sampler with the campaign
campaign.set_sampler(sampler)
Here, max_polynomial_order = 4 limits the maximum polynomial order of the individual interpolation stencils to 4.
Now we run all actions in the defined actions object, including running the code on the initial sampling plan:
###############################
# execute the defined actions #
###############################
campaign.execute().collate()
Now we set the following flags;
# number of samples used in the LEC check
n_mc_LEC = 50
# max number of LEC evaluations that can be run in parallel
max_LEC_jobs = 2
# max number of ENO stencils that can be computed in parallel
max_ENO_jobs = 2
# max number of functions evaluations that can be run in parallel
max_fjobs = 4
Here;
n_mc_LEC = 50: sets the number of times each local surrogate $w_j$ is evaluated in the LEC check to 50, with $j=1,\cdots,n_e$. Here $n_e$ is the number of simplex elements.max_LEC_jobs = 2: sets the number of LEC evaluations that are done in parallel to 4. Each simplex element can be given the LEC check independent of the others. The multiprocessing library is used for this.max_ENO_jobs = 2: sets the number of ENO stencils that are computed in parallel to 2. Again, the multiprocessing library is used for this.max_fjobs = 4: sets the number of new code evaluations per iteration of the SSC algorithm to 4. This means that the 4 simplex elements with the largest refinement metric $\bar{e}_j$ will get refined. Themax_f_jobscode evaluations can be done in parallel, and on supercomputers, using the VECMA/SEAVEA tools FabSim3 (see this recent article for a detailed overview of FabSim3) or QCG-PilotJob. In this tutorial the executing will be done sequentially on the localhost.
To refine the sampling plan, we need an SSCAnalysis object:
analysis = uq.analysis.SSCAnalysis(sampler = sampler, qoi_cols=['f'])
Here we update the surrogate (i.e. check the LEC condition and compute the ENO stencils), using the initial code samples (contained in the data_frame):
# get code samples
data_frame = campaign.get_collation_result()
# update surrogate by checking LEC condition and computing ENO stencils
analysis.update_surrogate('f', data_frame, n_mc_LEC=n_mc_LEC, max_LEC_jobs=max_LEC_jobs,
max_ENO_jobs=max_ENO_jobs)
# plot the sampling plan
analysis.plot_grid()
Loading samples... done Max. polynomial order allowed by n_s = 1 Checking LEC condition of 4 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1 1 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0
This shows the initial sampling plan, consisting of 4 simplices ($n_e = 4$) and 5 code samples ($n_s=5$). The code samples are located at the points of the simplices, Hence the initial sampling plan is created by placing one code sample in the middle, and one on each corner of the hypercube (in this case the square). Hence in the case of $n_\xi$ random inputs, the number of initial code samples is given by:
\begin{align*} n_s = 2^{n_\xi} + 1 \end{align*}The exponential increase the number of (initial) code samples is already an indication that the curse of dimensionality is present in the SSC method. As such, this method does not scale well to high $n_\xi$. Means to extend it to higher dimensions are discussed in reference [2] above.
For now, we stick with $n_\xi=2$. Note from the figure above that all elements are linear ($p_j=1$), since $n_s=5$ is not suffient to build higher-order interpolation stencils in 2D. Let us now execute our first refinement iteration:
def refine(campaign, analysis, n_refinements, max_fjobs, max_LEC_jobs = 4, n_mc_LEC = 50, max_ENO_jobs = 4):
"""
Refine the SSC sampling plan with max_fjobs new samples.
Parameters
----------
campaign : easyvvuq.campaign.Campaign
The EasyVVUQ campaign of the SSC sampler.
analysis : easyvvuq.analysis.ssc_analysis.SSCAnalysis
The SSC analysis object, used to refine the sampling plan.
n_refinements : int
The number of refinement iterations.
max_fjobs : int
The number of new code evaluations per refinement iteration.
max_LEC_jobs : int
The number of LEC checks to perform in parallel.
n_mc_LEC : int
The number of surrogate Monte Carlo samples used in each LEC check.
max_ENO_jobs : int
The number of ENO stencils to compute in parallel.
Returns
-------
None.
"""
for i in range(n_refinements):
# get code samples
data_frame = campaign.get_collation_result()
# update surrogate by checking LEC condition and computing ENO stencils
analysis.update_surrogate('f', data_frame, n_mc_LEC=n_mc_LEC, max_LEC_jobs=max_LEC_jobs,
max_ENO_jobs=max_ENO_jobs)
# refine the max_jobs simplices with the highest refinement measure
analysis.adapt_locally(max_fjobs)
# run the new code samples
campaign.execute().collate()
# do not forget to update the surrogate after the last code evaluations
data_frame = campaign.get_collation_result()
analysis.update_surrogate('f', data_frame)
# fix random seed
# np.random.seed(42)
# refine the sampling plan once
refine(campaign, analysis, 1, max_fjobs, n_mc_LEC)
analysis.plot_grid()
Loading samples... done Max. polynomial order allowed by n_s = 1 Checking LEC condition of 4 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1 1 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge refining edge refining edge Loading samples... done Max. polynomial order allowed by n_s = 2 Checking LEC condition of 8 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [2 2 2 1 1 1 1 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0
We can see that new code samples have been added, and that some simplices now have quadratic interpolation stencils, unless the LEC check reduced them back to 1. The refinement within the subsimplex is random, so if you rerun the analysis you will not get the exact same sampling plan, unless the random seed is fixed (uncomment the np.random.seed line above).
Let us now run several consecutive refinements:
refine(campaign, analysis, 15, max_fjobs, n_mc_LEC)
Loading samples... done Max. polynomial order allowed by n_s = 2 Checking LEC condition of 8 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1 1 1 1 1 1 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge refining edge Edge already refined, selecting another sample. Loading samples... done Max. polynomial order allowed by n_s = 3 Checking LEC condition of 16 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [2 2 1 1 1 1 1 1 1 1 1 1 2 2 1 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge Loading samples... done Max. polynomial order allowed by n_s = 4 Checking LEC condition of 23 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [2 2 1 1 1 2 2 2 1 2 1 1 2 1 1 1 2 2 2 1 1 2 2] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 30 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [2 2 2 2 2 1 2 2 2 1 2 2 2 1 1 1 1 1 1 1 1 2 2 2 2 1 1 1 1 2] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 36 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [3 3 2 3 2 3 3 3 2 2 3 3 3 3 3 1 2 1 1 1 1 2 2 2 2 1 1 1 2 2 2 1 1 3 2 2] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 41 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [3 3 2 3 2 3 3 3 2 1 3 3 3 3 1 1 1 2 2 2 1 1 1 1 2 2 2 2 2 1 1 1 1 2 1 1 2 2 2 3 2] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 49 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [3 3 3 3 3 3 3 3 3 3 3 3 1 1 3 3 3 2 1 3 3 3 2 2 2 1 1 3 3 3 3 3 1 1 1 3 2 2 2 1 1 1 1 1 2 1 1 2 2] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 57 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 1 4 3 3 4 4 3 1 4 1 3 3 4 4 3 1 2 2 1 2 2 2 2 1 1 2 1 2 1 1 1 3 1 3 3 3 3 3 4 4 4 4 4 4 2 2 2 4 4 1 1 1 1 3 3] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 65 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 4 3 4 4 1 1 3 1 1 3 3 3 4 4 2 2 2 2 3 1 4 4 4 3 3 3 3 3 4 4 4 3 3 1 1 1 3 3 1 3 3 1 2 2 2 4 4 2 2 2 2 1 2 1 2 1 1 3 1 1 1 1 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 72 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [3 3 3 3 3 3 1 3 1 1 4 4 4 1 3 2 3 3 2 1 3 3 3 4 4 4 4 4 3 3 3 4 4 3 3 1 4 4 1 1 2 2 2 3 3 3 3 1 1 2 2 3 3 3 3 3 1 4 1 1 4 4 4 1 1 1 1 1 2 2 1 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 78 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 3 4 4 4 4 4 4 3 1 1 4 4 4 2 3 3 4 1 4 4 4 4 4 4 4 4 4 4 4 4 1 4 4 1 1 2 2 2 4 4 3 3 1 1 1 3 3 3 3 3 3 4 4 3 3 4 4 4 1 1 1 4 1 1 4 4 4 1 1 1 1 3 3 1 1 1 3] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 86 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 3 4 4 4 4 4 4 3 1 1 3 4 4 2 4 4 3 3 4 1 1 4 4 4 4 4 4 4 4 4 4 4 3 3 1 2 2 1 3 3 3 1 1 1 3 3 3 3 4 4 4 4 4 4 4 1 3 1 1 3 1 4 3 3 1 1 1 2 1 1 1 1 1 4 4 4 3 1 1 1 1 1 4 4 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 94 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 2 4 4 4 4 4 4 2 4 3 4 4 4 4 4 2 4 4 3 2 2 1 1 4 1 4 4 4 1 4 4 4 4 4 4 4 4 3 1 2 2 1 3 3 3 3 3 1 1 3 1 1 1 2 1 1 1 4 3 3 3 3 1 1 1 3 1 1 3 4 1 1 1 3 3 4 4 4 4 1 1 1 3 3 1 1 1 3 4 1 1 4 3] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 102 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 4 4 4 4 4 4 4 4 3 4 4 4 4 4 2 4 4 3 4 4 4 1 1 4 1 4 4 4 1 4 4 4 4 4 4 4 4 3 1 2 2 1 3 3 3 3 1 1 3 3 1 1 2 3 2 1 3 4 3 3 1 1 1 1 3 3 3 1 1 4 4 1 1 1 3 3 4 4 4 4 3 2 1 1 1 1 3 3 3 3 1 1 3 1 4 1 4 3 1 1] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 109 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 4 4 4 4 4 4 4 4 4 4 4 4 2 4 4 3 4 4 1 1 4 4 4 4 1 4 4 4 4 4 4 4 4 1 2 2 1 3 3 3 3 1 4 1 3 1 3 1 1 4 4 4 3 4 4 1 1 3 4 3 3 3 3 1 1 1 1 3 4 4 1 1 1 4 4 4 4 4 4 4 1 1 1 1 3 1 1 1 1 1 3 3 4 4 4 1 4 1 4 2 2 4 1 1 2 3 3] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0 refining edge refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 115 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 2 4 4 1 1 4 4 4 4 1 4 4 4 4 4 4 4 4 1 3 2 4 4 1 4 1 3 1 4 1 1 4 4 4 1 2 3 4 4 2 1 3 1 2 3 2 4 4 4 3 3 1 1 4 4 4 2 2 1 1 4 4 4 4 4 4 2 2 4 1 2 4 1 1 3 4 1 1 1 2 1 4 4 4 4 4 4 1 4 2 4 2 2 4 1 1 2 3 3] Computing simplex probabilities... done, used 100000 samples. % of simplices with no samples = 0.0
analysis.plot_grid()
Where the discontinuity lies should now be clear through a relatively thin line of simplices with $p_j=1$ or $p_j=2$. Simplices with higher 4-th order interpolation stencils should appear away from the discontinutity, where the function is smooth.
We now plot the SSC surrogate;
plot_surrogate(analysis)
This should now be a much better approximation to the real function than the SC surrogate, free of the Runge phenomenom.
We can now compute the mean and variance with confidence. The SSC method uses simple Monte Carlo to compute these:
campaign.apply_analysis(analysis)
results = campaign.get_last_analysis()
results.describe()
Loading samples... done Computing mean and variance... done.
| f | |
|---|---|
| 0 | |
| mean | 0.896950 |
| std | 1.248267 |
| var | 1.558171 |
If you want to manually check the points and simplex elements, you can run the following
tri = sampler.get_Delaunay()
# indiviual points
print(tri.points)
# simplex indices
print(tri.simplices)
# simplex points
print(tri.points[tri.simplices])
[[0. 0. ] [0. 1. ] [1. 0. ] [1. 1. ] [0.5 0.5 ] [0.54064773 1. ] [0. 0.41521457] [1. 0.41644084] [0.4934277 0. ] [0.20274601 0.79725399] [0.75388736 0.75388736] [0.28215246 0.28742651] [0.60437879 0.24234286] [0.78459488 0.13194918] [0.38511521 0.85249404] [0.72686922 0.61002265] [1. 0.74412362] [0.13664381 0.47472514] [0.38019495 0.42948316] [0.81906525 0.4056478 ] [0. 0.68231194] [0.56551612 0.64490633] [1. 0.19311265] [0.19528432 0.61547185] [0.29793001 1. ] [0.32664753 0.72276204] [0.71601848 0. ] [0. 0.85285712] [0.17742196 1. ] [0.56502743 0.36974661] [0.87404743 0.21924869] [0.6395195 0.52793074] [0.76189818 0.24795324] [0.31876619 0.60327407] [0.11315931 0.86035105] [0.45924229 0.63530335] [0.37822795 0.55006621] [0.73339237 0.32661225] [0.68433183 0.43493004] [0.13067784 0.72687984] [0.87832922 0.34153458] [1. 0.31839895] [0.25939326 0.72319204] [0.06009197 0.79492669] [0.56985068 0.50294115] [0.2431041 0. ] [0.30835918 0.79230879] [0.10976166 0.68566957] [0. 0.17238511] [0.36370968 0.65979093] [0.66056549 0.37288498] [0.9300334 0.26332102] [0.60793171 0.42116799] [0.48552849 0.57878535] [0.45025676 0.56347211] [0.19270426 0.67340211] [0.26227833 0.63063678] [0.85621934 0.2776845 ] [0.73774888 0.39453742] [0.05619816 0.711593 ] [0.93788858 0.37820213] [0.84094386 0.52917209] [0. 0.76798322] [0.80857201 0.33809961] [0.37623423 0.60986567] [0.86398594 0. ] [1. 0.25064634] [0.16586081 0.72951474] [0.22186137 0.70675523]] [[11 17 6] [17 20 6] [48 11 6] [11 45 8] [45 48 0] [48 45 11] [12 11 8] [26 12 8] [17 47 20] [24 14 5] [ 5 10 3] [10 16 3] [36 18 4] [18 17 11] [12 13 32] [26 13 12] [27 34 1] [47 23 55] [23 47 17] [53 54 4] [54 36 4] [35 54 53] [33 23 17] [18 33 17] [33 18 36] [44 53 4] [21 31 15] [21 35 53] [44 21 53] [21 44 31] [10 21 15] [21 14 35] [14 21 5] [21 10 5] [22 65 2] [65 13 26] [34 28 1] [28 9 24] [ 9 28 34] [59 62 20] [47 59 20] [62 43 27] [43 34 27] [59 43 62] [14 25 35] [44 52 31] [52 44 4] [16 61 7] [61 10 15] [10 61 16] [66 51 22] [66 41 51] [60 41 7] [61 60 7] [60 61 19] [41 60 51] [30 65 22] [65 30 13] [51 30 22] [57 30 51] [13 30 32] [30 57 32] [46 14 24] [ 9 46 24] [46 25 14] [39 9 34] [43 39 34] [39 59 47] [39 43 59] [25 49 35] [23 56 55] [33 56 23] [49 56 33] [56 49 25] [56 68 55] [52 38 31] [31 38 15] [38 61 15] [18 29 4] [29 52 4] [12 29 11] [29 18 11] [40 57 51] [60 40 51] [40 60 19] [64 33 36] [64 49 33] [49 64 35] [64 54 35] [54 64 36] [39 67 9] [67 68 9] [68 67 55] [67 47 55] [67 39 47] [42 56 25] [42 68 56] [46 42 25] [42 46 9] [68 42 9] [61 58 19] [38 58 61] [58 63 19] [57 63 32] [63 40 19] [40 63 57] [29 50 52] [50 38 52] [50 58 38] [50 29 12] [63 37 32] [37 63 58] [50 37 58] [37 12 32] [37 50 12]] [[[0.28215246 0.28742651] [0.13664381 0.47472514] [0. 0.41521457]] [[0.13664381 0.47472514] [0. 0.68231194] [0. 0.41521457]] [[0. 0.17238511] [0.28215246 0.28742651] [0. 0.41521457]] [[0.28215246 0.28742651] [0.2431041 0. ] [0.4934277 0. ]] [[0.2431041 0. ] [0. 0.17238511] [0. 0. ]] [[0. 0.17238511] [0.2431041 0. ] [0.28215246 0.28742651]] [[0.60437879 0.24234286] [0.28215246 0.28742651] [0.4934277 0. ]] [[0.71601848 0. ] [0.60437879 0.24234286] [0.4934277 0. ]] [[0.13664381 0.47472514] [0.10976166 0.68566957] [0. 0.68231194]] [[0.29793001 1. ] [0.38511521 0.85249404] [0.54064773 1. ]] [[0.54064773 1. ] [0.75388736 0.75388736] [1. 1. ]] [[0.75388736 0.75388736] [1. 0.74412362] [1. 1. ]] [[0.37822795 0.55006621] [0.38019495 0.42948316] [0.5 0.5 ]] [[0.38019495 0.42948316] [0.13664381 0.47472514] [0.28215246 0.28742651]] [[0.60437879 0.24234286] [0.78459488 0.13194918] [0.76189818 0.24795324]] [[0.71601848 0. ] [0.78459488 0.13194918] [0.60437879 0.24234286]] [[0. 0.85285712] [0.11315931 0.86035105] [0. 1. ]] [[0.10976166 0.68566957] [0.19528432 0.61547185] [0.19270426 0.67340211]] [[0.19528432 0.61547185] [0.10976166 0.68566957] [0.13664381 0.47472514]] [[0.48552849 0.57878535] [0.45025676 0.56347211] [0.5 0.5 ]] [[0.45025676 0.56347211] [0.37822795 0.55006621] [0.5 0.5 ]] [[0.45924229 0.63530335] [0.45025676 0.56347211] [0.48552849 0.57878535]] [[0.31876619 0.60327407] [0.19528432 0.61547185] [0.13664381 0.47472514]] [[0.38019495 0.42948316] [0.31876619 0.60327407] [0.13664381 0.47472514]] [[0.31876619 0.60327407] [0.38019495 0.42948316] [0.37822795 0.55006621]] [[0.56985068 0.50294115] [0.48552849 0.57878535] [0.5 0.5 ]] [[0.56551612 0.64490633] [0.6395195 0.52793074] [0.72686922 0.61002265]] [[0.56551612 0.64490633] [0.45924229 0.63530335] [0.48552849 0.57878535]] [[0.56985068 0.50294115] [0.56551612 0.64490633] [0.48552849 0.57878535]] [[0.56551612 0.64490633] [0.56985068 0.50294115] [0.6395195 0.52793074]] [[0.75388736 0.75388736] [0.56551612 0.64490633] [0.72686922 0.61002265]] [[0.56551612 0.64490633] [0.38511521 0.85249404] [0.45924229 0.63530335]] [[0.38511521 0.85249404] [0.56551612 0.64490633] [0.54064773 1. ]] [[0.56551612 0.64490633] [0.75388736 0.75388736] [0.54064773 1. ]] [[1. 0.19311265] [0.86398594 0. ] [1. 0. ]] [[0.86398594 0. ] [0.78459488 0.13194918] [0.71601848 0. ]] [[0.11315931 0.86035105] [0.17742196 1. ] [0. 1. ]] [[0.17742196 1. ] [0.20274601 0.79725399] [0.29793001 1. ]] [[0.20274601 0.79725399] [0.17742196 1. ] [0.11315931 0.86035105]] [[0.05619816 0.711593 ] [0. 0.76798322] [0. 0.68231194]] [[0.10976166 0.68566957] [0.05619816 0.711593 ] [0. 0.68231194]] [[0. 0.76798322] [0.06009197 0.79492669] [0. 0.85285712]] [[0.06009197 0.79492669] [0.11315931 0.86035105] [0. 0.85285712]] [[0.05619816 0.711593 ] [0.06009197 0.79492669] [0. 0.76798322]] [[0.38511521 0.85249404] [0.32664753 0.72276204] [0.45924229 0.63530335]] [[0.56985068 0.50294115] [0.60793171 0.42116799] [0.6395195 0.52793074]] [[0.60793171 0.42116799] [0.56985068 0.50294115] [0.5 0.5 ]] [[1. 0.74412362] [0.84094386 0.52917209] [1. 0.41644084]] [[0.84094386 0.52917209] [0.75388736 0.75388736] [0.72686922 0.61002265]] [[0.75388736 0.75388736] [0.84094386 0.52917209] [1. 0.74412362]] [[1. 0.25064634] [0.9300334 0.26332102] [1. 0.19311265]] [[1. 0.25064634] [1. 0.31839895] [0.9300334 0.26332102]] [[0.93788858 0.37820213] [1. 0.31839895] [1. 0.41644084]] [[0.84094386 0.52917209] [0.93788858 0.37820213] [1. 0.41644084]] [[0.93788858 0.37820213] [0.84094386 0.52917209] [0.81906525 0.4056478 ]] [[1. 0.31839895] [0.93788858 0.37820213] [0.9300334 0.26332102]] [[0.87404743 0.21924869] [0.86398594 0. ] [1. 0.19311265]] [[0.86398594 0. ] [0.87404743 0.21924869] [0.78459488 0.13194918]] [[0.9300334 0.26332102] [0.87404743 0.21924869] [1. 0.19311265]] [[0.85621934 0.2776845 ] [0.87404743 0.21924869] [0.9300334 0.26332102]] [[0.78459488 0.13194918] [0.87404743 0.21924869] [0.76189818 0.24795324]] [[0.87404743 0.21924869] [0.85621934 0.2776845 ] [0.76189818 0.24795324]] [[0.30835918 0.79230879] [0.38511521 0.85249404] [0.29793001 1. ]] [[0.20274601 0.79725399] [0.30835918 0.79230879] [0.29793001 1. ]] [[0.30835918 0.79230879] [0.32664753 0.72276204] [0.38511521 0.85249404]] [[0.13067784 0.72687984] [0.20274601 0.79725399] [0.11315931 0.86035105]] [[0.06009197 0.79492669] [0.13067784 0.72687984] [0.11315931 0.86035105]] [[0.13067784 0.72687984] [0.05619816 0.711593 ] [0.10976166 0.68566957]] [[0.13067784 0.72687984] [0.06009197 0.79492669] [0.05619816 0.711593 ]] [[0.32664753 0.72276204] [0.36370968 0.65979093] [0.45924229 0.63530335]] [[0.19528432 0.61547185] [0.26227833 0.63063678] [0.19270426 0.67340211]] [[0.31876619 0.60327407] [0.26227833 0.63063678] [0.19528432 0.61547185]] [[0.36370968 0.65979093] [0.26227833 0.63063678] [0.31876619 0.60327407]] [[0.26227833 0.63063678] [0.36370968 0.65979093] [0.32664753 0.72276204]] [[0.26227833 0.63063678] [0.22186137 0.70675523] [0.19270426 0.67340211]] [[0.60793171 0.42116799] [0.68433183 0.43493004] [0.6395195 0.52793074]] [[0.6395195 0.52793074] [0.68433183 0.43493004] [0.72686922 0.61002265]] [[0.68433183 0.43493004] [0.84094386 0.52917209] [0.72686922 0.61002265]] [[0.38019495 0.42948316] [0.56502743 0.36974661] [0.5 0.5 ]] [[0.56502743 0.36974661] [0.60793171 0.42116799] [0.5 0.5 ]] [[0.60437879 0.24234286] [0.56502743 0.36974661] [0.28215246 0.28742651]] [[0.56502743 0.36974661] [0.38019495 0.42948316] [0.28215246 0.28742651]] [[0.87832922 0.34153458] [0.85621934 0.2776845 ] [0.9300334 0.26332102]] [[0.93788858 0.37820213] [0.87832922 0.34153458] [0.9300334 0.26332102]] [[0.87832922 0.34153458] [0.93788858 0.37820213] [0.81906525 0.4056478 ]] [[0.37623423 0.60986567] [0.31876619 0.60327407] [0.37822795 0.55006621]] [[0.37623423 0.60986567] [0.36370968 0.65979093] [0.31876619 0.60327407]] [[0.36370968 0.65979093] [0.37623423 0.60986567] [0.45924229 0.63530335]] [[0.37623423 0.60986567] [0.45025676 0.56347211] [0.45924229 0.63530335]] [[0.45025676 0.56347211] [0.37623423 0.60986567] [0.37822795 0.55006621]] [[0.13067784 0.72687984] [0.16586081 0.72951474] [0.20274601 0.79725399]] [[0.16586081 0.72951474] [0.22186137 0.70675523] [0.20274601 0.79725399]] [[0.22186137 0.70675523] [0.16586081 0.72951474] [0.19270426 0.67340211]] [[0.16586081 0.72951474] [0.10976166 0.68566957] [0.19270426 0.67340211]] [[0.16586081 0.72951474] [0.13067784 0.72687984] [0.10976166 0.68566957]] [[0.25939326 0.72319204] [0.26227833 0.63063678] [0.32664753 0.72276204]] [[0.25939326 0.72319204] [0.22186137 0.70675523] [0.26227833 0.63063678]] [[0.30835918 0.79230879] [0.25939326 0.72319204] [0.32664753 0.72276204]] [[0.25939326 0.72319204] [0.30835918 0.79230879] [0.20274601 0.79725399]] [[0.22186137 0.70675523] [0.25939326 0.72319204] [0.20274601 0.79725399]] [[0.84094386 0.52917209] [0.73774888 0.39453742] [0.81906525 0.4056478 ]] [[0.68433183 0.43493004] [0.73774888 0.39453742] [0.84094386 0.52917209]] [[0.73774888 0.39453742] [0.80857201 0.33809961] [0.81906525 0.4056478 ]] [[0.85621934 0.2776845 ] [0.80857201 0.33809961] [0.76189818 0.24795324]] [[0.80857201 0.33809961] [0.87832922 0.34153458] [0.81906525 0.4056478 ]] [[0.87832922 0.34153458] [0.80857201 0.33809961] [0.85621934 0.2776845 ]] [[0.56502743 0.36974661] [0.66056549 0.37288498] [0.60793171 0.42116799]] [[0.66056549 0.37288498] [0.68433183 0.43493004] [0.60793171 0.42116799]] [[0.66056549 0.37288498] [0.73774888 0.39453742] [0.68433183 0.43493004]] [[0.66056549 0.37288498] [0.56502743 0.36974661] [0.60437879 0.24234286]] [[0.80857201 0.33809961] [0.73339237 0.32661225] [0.76189818 0.24795324]] [[0.73339237 0.32661225] [0.80857201 0.33809961] [0.73774888 0.39453742]] [[0.66056549 0.37288498] [0.73339237 0.32661225] [0.73774888 0.39453742]] [[0.73339237 0.32661225] [0.60437879 0.24234286] [0.76189818 0.24795324]] [[0.73339237 0.32661225] [0.66056549 0.37288498] [0.60437879 0.24234286]]]
One dimensional problems¶
For 2 or more random inputs the SSCSampler uses the SciPy Delaunay triangulation. This does not work if there is only 1 random input. You can still use the SSCSampler though, as under the hood the Delaunay triangulation (sampler.tri) is replaced by a 1D grid that mimics the SciPy Delaunay properties and subroutines that are expected by the SSCSampler.
This is illustrated below, where we remove $\xi_2$ as a random variable. The rest is unchanged.
campaign = get_campaign(n_xi = 1)
# redefine vary as 1D, this fixes xi_2 to its default value defined in params
vary = {'xi_1' : cp.Uniform(0, 1)}
##################
# Select sampler #
##################
sampler = uq.sampling.SSCSampler(vary=vary, max_polynomial_order=4)
# Associate the sampler with the campaign
campaign.set_sampler(sampler)
###############################
# execute the defined actions #
###############################
campaign.execute().collate()
analysis = uq.analysis.SSCAnalysis(sampler = sampler, qoi_cols=['f'])
refine(campaign, analysis, 5, max_fjobs, n_mc_LEC)
Loading samples... done Max. polynomial order allowed by n_s = 2 Checking LEC condition of 2 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1] Computing simplex probabilities... done, used 10000 samples. % of simplices with no samples = 0.0 refining edge refining edge Loading samples... done Max. polynomial order allowed by n_s = 4 Checking LEC condition of 4 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [2 2 1 1] Computing simplex probabilities... done, used 10000 samples. % of simplices with no samples = 0.0 refining edge refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 8 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 4 4 4 4 4 4] Computing simplex probabilities... done, used 10000 samples. % of simplices with no samples = 0.0 Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 12 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 4 4 4 4 3 4 4 4 4 4] Computing simplex probabilities... done, used 10000 samples. % of simplices with no samples = 0.0 refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 16 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 4 4 4 4 4 4 4 1 3 4 4 4 4 4] Computing simplex probabilities... done, used 10000 samples. % of simplices with no samples = 0.0 refining edge Loading samples... done Max. polynomial order set by hand to = 4 Checking LEC condition of 20 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [4 4 4 4 4 4 4 4 4 2 1 4 4 4 4 4 4 4 4 4] Computing simplex probabilities... done, used 10000 samples. % of simplices with no samples = 0.0
A plot of the 1D grid and the polynomial order of each element. We can already see where the LEC condition lowered the polynomial order to avoid oscillations at the discontinuity.
analysis.plot_grid()
The 1D surrogate:
%matplotlib inline
plot_surrogate(analysis)
Higher dimensional problems¶
As mentioned, the curse of dimensionality is present in the SSC method. Not only due to the exponential increase in the number of code evaluations, but also in the cost of building the surrogate [2]. This is related to the exponential increase in the number of simplex elements in the Delaunay triangulation. If you are running this one a node with a large number of cores, you can increase max_LEC_jobs and max_ENO_jobs to get some speedup.
This can be seen by executing the 5D example below. You will note that checking the LEC condition and computing the ENO stencils takes considerably longer. Methods for extending the SSC sampler to dimensions higher than 5 are discussed in reference [2] above (e.g. HDMR expansions). These will be implemented in a later stage.
n_xi = 5
campaign = get_campaign(n_xi = n_xi)
# redefine vary
for i in range(n_xi):
vary['xi_%d' % (i + 1)] = cp.Uniform(0, 1)
##################
# Select sampler #
##################
sampler = uq.sampling.SSCSampler(vary=vary, max_polynomial_order=4)
# Associate the sampler with the campaign
campaign.set_sampler(sampler)
###############################
# execute the defined actions #
###############################
campaign.execute().collate()
analysis = uq.analysis.SSCAnalysis(sampler = sampler, qoi_cols=['f'])
max_fjobs = 5
# increase these a little to speed things up if you have available cores
max_LEC_jobs = 8
max_ENO_jobs = 8
refine(campaign, analysis, 5, max_fjobs, max_LEC_jobs=max_LEC_jobs, max_ENO_jobs=max_ENO_jobs)
Loading samples... done Max. polynomial order allowed by n_s = 2 Checking LEC condition of 210 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 1 1 1 1 1 1 1 1 1 1 1 1 2 1 2 1 2 2 2 2 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2] Computing simplex probabilities... done, used 10000000 samples. % of simplices with no samples = 0.0 refining edge refining edge refining edge refining edge refining edge Loading samples... done Max. polynomial order allowed by n_s = 2 Checking LEC condition of 418 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 1 1 1 1 2 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 1 1 2 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 2 1 1 1 1 1 1 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 1 1 1 1 2 1 1 1 1 1 2 1 2 1 2 1 2 1 1 1 1 1 2 1 1 2 1 1 1 1 1 1 2 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 1 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1 2 1 1 2 2 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1] Computing simplex probabilities... done, used 10000000 samples. % of simplices with no samples = 0.0 refining edge refining edge refining edge refining edge refining edge Edge already refined, selecting another sample. refining edge Loading samples... done Max. polynomial order allowed by n_s = 2 Checking LEC condition of 535 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1] Computing simplex probabilities... done, used 10000000 samples. % of simplices with no samples = 0.0 refining edge refining edge refining edge refining edge refining edge Loading samples... done Max. polynomial order allowed by n_s = 2 Checking LEC condition of 694 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1] Computing simplex probabilities... done, used 10000000 samples. % of simplices with no samples = 0.0 refining edge refining edge refining edge refining edge refining edge Edge already refined, selecting another sample. refining edge Loading samples... done Max. polynomial order allowed by n_s = 2 Checking LEC condition of 916 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1] Computing simplex probabilities... done, used 10000000 samples. % of simplices with no samples = 0.0 refining edge refining edge Edge already refined, selecting another sample. refining edge refining edge Edge already refined, selecting another sample. refining edge refining edge refining edge Loading samples... done Max. polynomial order allowed by n_s = 3 Checking LEC condition of 1044 stencils... done. Computing ENO stencils... done. Updated polynomials orders = [1 1 1 ... 1 1 1] Computing simplex probabilities... done, used 10000000 samples. % of simplices with no samples = 0.0