Polynomial Chaos Surrogate

This example explains how to train and utilize a PolynomialChaos surrogate model. For detailed information on the polynomial chaos method including the mathematics and implementation see PolynomialChaos and QuadratureSampler. For general information on training and evaluating a surrogate model see Training a surrogate model and Evaluating a surrogate model.

Overview

Polynomial chaos (PC) is a stochastic element method. Its formulation is similar to finite element method, but instead of discretizing spatial dimensions, PC discretizes random variables, i.e. uncertain parameters. Within the training phase, PC creates a functional of a quantity of interest (QoI) that is dependent on the random variables. This function consists of a finite sum of orthogonal polynomials with coefficients, these coefficients are what the training phase is attempting to compute. Another difference from finite elements is that these polynomials do not necessarily have compact support and increasing the terms in the sum increases the polynomial order. Please see PolynomialChaos for more details of the polynomial chaos method.

Model Problem

This example uses a one-dimensional heat conduction problem as the full-order model which has certain uncertain parameters. The model equation is as follows:

The quantities of interest are the average and maximum temperature:

Parameter Uncertainty

For demonstration, each of these parameters will have two types of probability distributions: Uniform () and Normal (). Where and are the max and minimum bounds of the uniform distribution, respectively. And and are the mean and standard deviation of the normal distribution, respectively.

The uncertain parameters for this model problem are:

ParameterSymbolUniformNormal
Conductivity
Volumetric Heat Source
Domain Size
Right Boundary Temperature

Analytical Solutions

This simple model problem has analytical descriptions for the field temperature, average temperature, and maximum temperature:

With the quadratic feature of the field temperature, using quadratic elements in the discretization will actually yield the exact solution.

Input File

Below is the input file used to solve the one-dimensional heat conduction model.

[Mesh<<<{"href": "../../../syntax/Mesh/index.html"}>>>]
  type = GeneratedMesh
  dim = 1
  nx = 100
  xmax = 1
  elem_type = EDGE3
[]

[Variables<<<{"href": "../../../syntax/Variables/index.html"}>>>]
  [T]
    order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = SECOND
    family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = LAGRANGE
  []
[]

[Kernels<<<{"href": "../../../syntax/Kernels/index.html"}>>>]
  [diffusion]
    type = MatDiffusion<<<{"description": "Diffusion equation Kernel that takes an isotropic Diffusivity from a material property", "href": "../../../source/kernels/MatDiffusion.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    diffusivity<<<{"description": "The diffusivity value or material property"}>>> = k
  []
  [source]
    type = BodyForce<<<{"description": "Demonstrates the multiple ways that scalar values can be introduced into kernels, e.g. (controllable) constants, functions, and postprocessors. Implements the weak form $(\\psi_i, -f)$.", "href": "../../../source/kernels/BodyForce.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    value<<<{"description": "Coefficient to multiply by the body force term"}>>> = 1.0
  []
[]

[Materials<<<{"href": "../../../syntax/Materials/index.html"}>>>]
  [conductivity]
    type = GenericConstantMaterial<<<{"description": "Declares material properties based on names and values prescribed by input parameters.", "href": "../../../source/materials/GenericConstantMaterial.html"}>>>
    prop_names<<<{"description": "The names of the properties this material will have"}>>> = k
    prop_values<<<{"description": "The values associated with the named properties"}>>> = 2.0
  []
[]

[BCs<<<{"href": "../../../syntax/BCs/index.html"}>>>]
  [right]
    type = DirichletBC<<<{"description": "Imposes the essential boundary condition $u=g$, where $g$ is a constant, controllable value.", "href": "../../../source/bcs/DirichletBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = right
    value<<<{"description": "Value of the BC"}>>> = 300
  []
[]

[Executioner<<<{"href": "../../../syntax/Executioner/index.html"}>>>]
  type = Steady
  solve_type = PJFNK
  petsc_options_iname = '-pc_type -pc_hypre_type'
  petsc_options_value = 'hypre boomeramg'
[]

[Postprocessors<<<{"href": "../../../syntax/Postprocessors/index.html"}>>>]
  [avg]
    type = AverageNodalVariableValue<<<{"description": "Computes the average value of a field by sampling all nodal solutions on the domain or within a subdomain", "href": "../../../source/postprocessors/AverageNodalVariableValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on"}>>> = T
  []
  [max]
    type = NodalExtremeValue<<<{"description": "Finds either the min or max elemental value of a variable over the domain.", "href": "../../../source/postprocessors/NodalExtremeValue.html"}>>>
    variable<<<{"description": "The name of the variable that this postprocessor operates on"}>>> = T
    value_type<<<{"description": "Type of extreme value to return. 'max' returns the maximum value. 'min' returns the minimum value. 'max_abs' returns the maximum of the absolute value."}>>> = max
  []
[]

[Outputs<<<{"href": "../../../syntax/Mastodon/Outputs/index.html"}>>>]
[]
(moose/modules/stochastic_tools/examples/surrogates/sub.i)

With this input the uncertain parameters are defined as:

  1. Materials/conductivity/prop_values

  2. Kernels/source/value

  3. Mesh/xmax

  4. BCs/right/value

These values in the sub.i file are arbitrary since the stochastic master app will be modifying them.

Training

This section describes how to set up an input for a PolynomialChaosTrainer. There are three aspects that need to be known before setting up the input file.

  1. A sub app needs to be built, defining the physics of the problem, this example uses a heat conduction model described in the previous section.

  2. Uncertain parameters need to be identified with a defined probability distribution.

  3. Quantities of interest need to be defined, for now, these come in the form of postprocessors from the sub app.

Omitting Solve

Any input file in MOOSE needs to include a Mesh, Variables, and Executioner block. However, the stochastic master app does not actually create or solve a system. So the StochasticToolsAction builds a minimal model to satisfy these requirements:

[StochasticTools<<<{"href": "../../../syntax/StochasticTools/index.html"}>>>]
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform_mc.i)

Parameter Uncertainty Distributions

A required aspect of the polynomial chaos method is a definition of the probability distribution of the uncertain parameters. These distributions define the type of polynomial used for the PC expansion. This example shows the two types of distributions available for PC training:

Uniform distribution for each parameter

[Distributions<<<{"href": "../../../syntax/Distributions/index.html"}>>>]
  [k_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 1
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 10
  []
  [q_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 9000
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 11000
  []
  [L_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 0.01
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 0.05
  []
  [Tinf_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 290
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 310
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform_mc.i)

Normal distribution for each parameter

[Distributions<<<{"href": "../../../syntax/Distributions/index.html"}>>>]
  [k_dist]
    type = Normal<<<{"description": "Normal distribution", "href": "../../../source/distributions/Normal.html"}>>>
    mean<<<{"description": "Mean (or expectation) of the distribution."}>>> = 5
    standard_deviation<<<{"description": "Standard deviation of the distribution "}>>> = 2
  []
  [q_dist]
    type = Normal<<<{"description": "Normal distribution", "href": "../../../source/distributions/Normal.html"}>>>
    mean<<<{"description": "Mean (or expectation) of the distribution."}>>> = 10000
    standard_deviation<<<{"description": "Standard deviation of the distribution "}>>> = 500
  []
  [L_dist]
    type = Normal<<<{"description": "Normal distribution", "href": "../../../source/distributions/Normal.html"}>>>
    mean<<<{"description": "Mean (or expectation) of the distribution."}>>> = 0.03
    standard_deviation<<<{"description": "Standard deviation of the distribution "}>>> = 0.01
  []
  [Tinf_dist]
    type = Normal<<<{"description": "Normal distribution", "href": "../../../source/distributions/Normal.html"}>>>
    mean<<<{"description": "Mean (or expectation) of the distribution."}>>> = 300
    standard_deviation<<<{"description": "Standard deviation of the distribution "}>>> = 10
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_normal_mc.i)

Training Sampler

Non-intrusive polynomial chaos training relies on sampling the full-order model (i.e. sub app) at a number of sample points. There are two general types of sampling techniques for PC: Monte Carlo and numerical quadrature.

Monte Carlo Sampling

Monte Carlo sampling is the most intuitive type of sampling and can be done using the MonteCarlo or LatinHypercube samplers in the stochastic tools module.

[Samplers<<<{"href": "../../../syntax/Samplers/index.html"}>>>]
  [sample]
    type = MonteCarlo<<<{"description": "Monte Carlo Sampler.", "href": "../../../source/samplers/MonteCarloSampler.html"}>>>
    num_rows<<<{"description": "The number of rows per matrix to generate."}>>> = 10000
    distributions<<<{"description": "The distribution names to be sampled, the number of distributions provided defines the number of columns per matrix."}>>> = 'k_dist q_dist L_dist Tinf_dist'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = PRE_MULTIAPP_SETUP
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform_mc.i)

Quadrature Sampling

Because PC expansion uses very specific types of polynomials, using numerical quadrature is typically a more precise approach to training a PC model. The QuadratureSampler implements numerical quadrature for uniform and normal distributions. There are three different techniques for building a multidimensional quadrature grid: Cartesian grid, Smolyak sparse grid, and Clenshaw-Curtis sparse grid. These grids are specified by the "sparse_grid" parameter: none, smolyak, and clenshaw-curtis. Using sparse grids can significantly reduce the number of training points while retaining accuracy, see QuadratureSampler for more details. The number of sample points is ultimately defined by the "order" parameter. It is generally advisable that this parameter be set to the same value as the "order" parameter in PolynomialChaosTrainer.

[Samplers<<<{"href": "../../../syntax/Samplers/index.html"}>>>]
  [sample]
    type = Quadrature<<<{"description": "Quadrature sampler for Polynomial Chaos.", "href": "../../../source/samplers/QuadratureSampler.html"}>>>
    order<<<{"description": "Specify the maximum order of the polynomials in the expansion."}>>> = 10
    distributions<<<{"description": "The distribution names to be sampled, the number of distributions provided defines the number of columns per matrix and their type defines the quadrature."}>>> = 'k_dist q_dist L_dist Tinf_dist'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = PRE_MULTIAPP_SETUP
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform_quad.i)

Running Sub Application

Running the sub app and transferring data back and forth for PC is exactly the same as any other type of surrogate. See Training a surrogate model on more details on setting up this part of the input file.

[MultiApps<<<{"href": "../../../syntax/MultiApps/index.html"}>>>]
  [sub]
    type = SamplerFullSolveMultiApp<<<{"description": "Creates a full-solve type sub-application for each row of each Sampler matrix.", "href": "../../../source/multiapps/SamplerFullSolveMultiApp.html"}>>>
    input_files<<<{"description": "The input file for each App.  If this parameter only contains one input file it will be used for all of the Apps.  When using 'positions_from_file' it is also admissable to provide one input_file per file."}>>> = sub.i
    sampler<<<{"description": "The Sampler object to utilize for creating MultiApps."}>>> = sample
  []
[]

[Controls<<<{"href": "../../../syntax/Controls/index.html"}>>>]
  [cmdline]
    type = MultiAppSamplerControl<<<{"description": "Control for modifying the command line arguments of MultiApps.", "href": "../../../source/controls/MultiAppSamplerControl.html"}>>>
    multi_app<<<{"description": "The name of the MultiApp to control."}>>> = sub
    sampler<<<{"description": "The Sampler object to utilize for altering the command line options of the MultiApp."}>>> = sample
    param_names<<<{"description": "The names of the command line parameters to set via the sampled data."}>>> = 'Materials/conductivity/prop_values Kernels/source/value Mesh/xmax BCs/right/value'
  []
[]

[Transfers<<<{"href": "../../../syntax/Transfers/index.html"}>>>]
  [data]
    type = SamplerReporterTransfer<<<{"description": "Transfers data from Reporters on the sub-application to a StochasticReporter on the main application.", "href": "../../../source/transfers/SamplerReporterTransfer.html"}>>>
    from_multi_app<<<{"description": "The name of the MultiApp to receive data from"}>>> = sub
    sampler<<<{"description": "A the Sampler object that Transfer is associated.."}>>> = sample
    stochastic_reporter<<<{"description": "The name of the StochasticReporter object to transfer values to."}>>> = results
    from_reporter<<<{"description": "The name(s) of the Reporter(s) on the sub-app to transfer from."}>>> = 'avg/value max/value'
  []
[]

[Reporters<<<{"href": "../../../syntax/Reporters/index.html"}>>>]
  [results]
    type = StochasticReporter<<<{"description": "Storage container for stochastic simulation results coming from Reporters.", "href": "../../../source/reporters/StochasticReporter.html"}>>>
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform_quad.i)

Trainers

The PC training object is defined within the Trainers block. The required parameters for a PC trainer are:

  • "distributions" specify the type of polynomials to use for the expansion, it is very import that these distributions match the distributions given to the sampler.

  • "sampler" is the object that will provide sample points which are given to the sub app during execution.

  • "response" specifies the result vector for storing the computed values.

  • "order" defines the maximum order of the PC expansion, this parameter ultimately defines the accuracy and complexity of the surrogate model.

[Trainers<<<{"href": "../../../syntax/Trainers/index.html"}>>>]
  [poly_chaos_avg]
    type = PolynomialChaosTrainer<<<{"description": "Computes and evaluates polynomial chaos surrogate model.", "href": "../../../source/trainers/PolynomialChaosTrainer.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    order<<<{"description": "Maximum polynomial order."}>>> = 10
    regression_type<<<{"description": "The type of regression to perform for finding polynomial coefficents."}>>> = integration
    distributions<<<{"description": "Names of the distributions samples were taken from."}>>> = 'k_dist q_dist L_dist Tinf_dist'
    sampler<<<{"description": "Sampler used to create predictor and response data."}>>> = sample
    response<<<{"description": "Reporter value of response results, can be vpp with <vpp_name>/<vector_name> or sampler column with 'sampler/col_<index>'."}>>> = results/data:avg:value
  []
  [poly_chaos_max]
    type = PolynomialChaosTrainer<<<{"description": "Computes and evaluates polynomial chaos surrogate model.", "href": "../../../source/trainers/PolynomialChaosTrainer.html"}>>>
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = timestep_end
    order<<<{"description": "Maximum polynomial order."}>>> = 10
    regression_type<<<{"description": "The type of regression to perform for finding polynomial coefficents."}>>> = integration
    distributions<<<{"description": "Names of the distributions samples were taken from."}>>> = 'k_dist q_dist L_dist Tinf_dist'
    sampler<<<{"description": "Sampler used to create predictor and response data."}>>> = sample
    response<<<{"description": "Reporter value of response results, can be vpp with <vpp_name>/<vector_name> or sampler column with 'sampler/col_<index>'."}>>> = results/data:max:value
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform_mc.i)

Outputting Training Data

Outputting the data after training a PC trainer is exactly the same as outputting data for any other type of trainer. See Training a surrogate model on more details on setting up this part of the input file.

[Outputs<<<{"href": "../../../syntax/Mastodon/Outputs/index.html"}>>>]
  file_base<<<{"description": "Common file base name to be utilized with all output objects"}>>> = poly_chaos_training
  [out]
    type = SurrogateTrainerOutput<<<{"description": "Output for trained surrogate model data.", "href": "../../../source/outputs/SurrogateTrainerOutput.html"}>>>
    trainers<<<{"description": "A list of SurrogateTrainer objects to output."}>>> = 'poly_chaos_avg poly_chaos_max'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = FINAL
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform_mc.i)

Example Input Files

Evaluation and Statistics

This section will go over how to set up an input file to evaluate and perform statistical analysis on a trained polynomial chaos surrogate. The polynomial chaos surrogate is unique in that it is able to compute statistical moments and sensitivities analytically. Specific postprocessors are available to compute these quantities without the need for sampling.

Omitting Solve

Any input file in MOOSE needs to include a Mesh, Variables, and Executioner block. However, the stochastic master app does not actually create or solve a system. So the StochasticToolsAction builds a minimal model to satisfy these requirements:

[StochasticTools<<<{"href": "../../../syntax/StochasticTools/index.html"}>>>]
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform.i)

Evaluation

Evaluating a polynomial chaos surrogate model is exactly the same as any other type of surrogate, see Evaluating a surrogate model for more details.

Loading training data

[Surrogates<<<{"href": "../../../syntax/Surrogates/index.html"}>>>]
  [poly_chaos_avg]
    type = PolynomialChaos<<<{"description": "Computes and evaluates polynomial chaos surrogate model.", "href": "../../../source/surrogates/PolynomialChaos.html"}>>>
    filename<<<{"description": "The name of the file which will be associated with the saved/loaded data."}>>> = 'poly_chaos_training_poly_chaos_avg.rd'
  []
  [poly_chaos_max]
    type = PolynomialChaos<<<{"description": "Computes and evaluates polynomial chaos surrogate model.", "href": "../../../source/surrogates/PolynomialChaos.html"}>>>
    filename<<<{"description": "The name of the file which will be associated with the saved/loaded data."}>>> = 'poly_chaos_training_poly_chaos_max.rd'
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform.i)

Defining sampler for evaluation – Uniform distributions

[Distributions<<<{"href": "../../../syntax/Distributions/index.html"}>>>]
  [k_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 1
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 10
  []
  [q_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 9000
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 11000
  []
  [L_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 0.01
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 0.05
  []
  [Tinf_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 290
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 310
  []
[]

[Samplers<<<{"href": "../../../syntax/Samplers/index.html"}>>>]
  [sample]
    type = MonteCarlo<<<{"description": "Monte Carlo Sampler.", "href": "../../../source/samplers/MonteCarloSampler.html"}>>>
    num_rows<<<{"description": "The number of rows per matrix to generate."}>>> = 100000
    distributions<<<{"description": "The distribution names to be sampled, the number of distributions provided defines the number of columns per matrix."}>>> = 'k_dist q_dist L_dist Tinf_dist'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = initial
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform.i)

Defining sampler for evaluation – Normal distributions

[Distributions<<<{"href": "../../../syntax/Distributions/index.html"}>>>]
  [k_dist]
    type = Normal<<<{"description": "Normal distribution", "href": "../../../source/distributions/Normal.html"}>>>
    mean<<<{"description": "Mean (or expectation) of the distribution."}>>> = 5
    standard_deviation<<<{"description": "Standard deviation of the distribution "}>>> = 2
  []
  [q_dist]
    type = Normal<<<{"description": "Normal distribution", "href": "../../../source/distributions/Normal.html"}>>>
    mean<<<{"description": "Mean (or expectation) of the distribution."}>>> = 10000
    standard_deviation<<<{"description": "Standard deviation of the distribution "}>>> = 500
  []
  [L_dist]
    type = Normal<<<{"description": "Normal distribution", "href": "../../../source/distributions/Normal.html"}>>>
    mean<<<{"description": "Mean (or expectation) of the distribution."}>>> = 0.03
    standard_deviation<<<{"description": "Standard deviation of the distribution "}>>> = 0.01
  []
  [Tinf_dist]
    type = Normal<<<{"description": "Normal distribution", "href": "../../../source/distributions/Normal.html"}>>>
    mean<<<{"description": "Mean (or expectation) of the distribution."}>>> = 300
    standard_deviation<<<{"description": "Standard deviation of the distribution "}>>> = 10
  []
[]

[Samplers<<<{"href": "../../../syntax/Samplers/index.html"}>>>]
  [sample]
    type = MonteCarlo<<<{"description": "Monte Carlo Sampler.", "href": "../../../source/samplers/MonteCarloSampler.html"}>>>
    num_rows<<<{"description": "The number of rows per matrix to generate."}>>> = 100000
    distributions<<<{"description": "The distribution names to be sampled, the number of distributions provided defines the number of columns per matrix."}>>> = 'k_dist q_dist L_dist Tinf_dist'
    execute_on<<<{"description": "The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html."}>>> = initial
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_normal.i)

Evaluating surrogate with EvaluateSurrogate

[Reporters<<<{"href": "../../../syntax/Reporters/index.html"}>>>]
  [samp]
    type = EvaluateSurrogate<<<{"description": "Tool for sampling surrogate models.", "href": "../../../source/reporters/EvaluateSurrogate.html"}>>>
    model<<<{"description": "Name of surrogate models."}>>> = 'poly_chaos_avg poly_chaos_max'
    sampler<<<{"description": "Sampler to use for evaluating surrogate models."}>>> = sample
    parallel_type<<<{"description": "This parameter will determine how the stochastic data is gathered. It is common for outputting purposes that this parameter be set to ROOT, otherwise, many files will be produced showing the values on each processor. However, if there are lot of samples, gathering on root may be memory restrictive."}>>> = ROOT
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform.i)

Statistics and Sensitivities

PolynomialChaos has the ability to compute statistical moments, local sensitivity, and sobol indices analytically. These are computed using PolynomialChaosReporter. The implemented moments include mean, standard deviation, skewness, and kurtosis. The type of moment is specified by the "statistics" parameter, which can be a combination of mean, stddev, skewness, and kurtosis. This example computes mean and standard deviation. It should be noted that computing skewness and kurtosis can be very computationally demanding using the analytical technique. Therefore, it might be better to sample the surrogate model and compute these quantities with the result.

The points to compute the sensitivity can be defined explicitly by specifying the "local_sensitivity_points" parameter, and/or using a sampler by specifying the "local_sensitivity_sampler". It is vital that the columns of the points match with the distributions defined from the original trainer.

Sobol statistics are a metric of the global sensitivity of each parameter, or a combination of parameters, see SobolReporter for more details. These are computed by setting the "include_sobol" parameter to true. This will compute total, first-, and second-order Sobol indices.

[Reporters<<<{"href": "../../../syntax/Reporters/index.html"}>>>]
  [stats]
    type = PolynomialChaosReporter<<<{"description": "Tool for extracting data from PolynomialChaos surrogates and computing statistics.", "href": "../../../source/reporters/PolynomialChaosReporter.html"}>>>
    pc_name<<<{"description": "Name(s) of PolynomialChaos surrogate object(s)."}>>> = 'poly_chaos_avg poly_chaos_max'
    statistics<<<{"description": "Statistics to compute."}>>> = 'mean stddev'
    local_sensitivity_points<<<{"description": "Points for each polynomial chaos surrogate specifying desired location of sensitivity measurement."}>>> = '5 10000 0.03 300; 5 10000 0.03 300'
    include_sobol<<<{"description": "True to compute Sobol indices."}>>> = true
  []
[]
(moose/modules/stochastic_tools/examples/surrogates/poly_chaos_uniform.i)

Example Input Files

Results and Analysis

In this section, the results of training and evaluation of the surrogate model as described in the previous sections is shown. A convergence study of the polynomial chaos training, including analysis of sampling and polynomial order, is also provided.

Evaluation Results

Figure 1 and Figure 2 show the resulting probability distributions from sampling each surrogate with 100,000 points.

Figure 1: Temperature distributions with uniform parameter distribution

Figure 2: Temperature distribution with normal parameter distribution

Local Sensitivity Results

Table 1 shows the results of computing the local sensitivity of each parameter for every surrogate. For reference, the sensitivity is defined as

where is either or , is , , , or , and the point tested () is

Table 1: Local sensitivity results

QoIDistribution Sensitivity Sensitivity Sensitivity Sensitivity
Uniform-0.0019440.0019850.0039760.9977
Normal-0.0020090.0016670.0033030.9983
Uniform-0.0029180.0029790.0059660.9967
Normal-0.0028970.0027360.0054370.9973

Sobol Statistics Results

Figure 3 shows the results of computing Sobol statistics for average temperature and uniform parameter distributions (the other surrogates produced very similar statistics). In the plot, the x axis represents the first index subscript and the y axis represents the second. As can be seen, the statistics are symmetric, i.e. it does not matter the order of the index subscripts.

Figure 3: Average temperature Sobol statistics results uniform parameter distribution

Sampling Analysis

Here the effect of various sampling techniques on the training of a polynomial chaos surrogate is shown. Three different types of samplers were used: Monte Carlo, tensor quadrature, and Smolyak sparse quadrature. Monte Carlo was run with increasing number of samples and the quadrature was run with increasing order. Note, the polynomial expansion order was set constant. The performance metric is the error in the mean and standard deviation versus the number of total training points. Figure 4-Figure 7 show the results of this study for each quantity of interest and each parameter distribution. The uniform distributions (Figure 4 and Figure 5) have an analytical mean and standard deviation, which is why there is nice convergence. However, the normal distributions (Figure 6 and Figure 7) use numerical integration to compute the reference mean and standard deviation, which is why there is a more spurious convergence at low errors.

Figure 4: Sample convergence for moments of average temperature with uniform parameter distributions

Figure 5: Sample convergence for moments of maximum temperature with uniform parameter distributions

Figure 6: Sample convergence for moments of average temperature with normal parameter distributions

Figure 7: Sample convergence for moments of maximum temperature with normal parameter distributions

Polynomial Order Analysis

The effect of polynomial expansion order on the accuracy of the surrogate model is shown here. Again, the performance metric is the error in the mean and standard deviation versus the number of total training points. The maximum polynomial order is ranged from 3 to 9, where the tensor quadrature sampler order matches (which is generally advisable in any case). Figure 8 shows the result of the convergence study with a uniform distribution. Normal distribution was omitted since there is no analytical reference.

Figure 8: Polynomial order convergence for moments of maximum and average temperature with uniform parameter distributions