Polynomial Regression Surrogate

This example is meant to demonstrate how a polynomial regression based surrogate model is trained and used on a parametric problem. Additionally, the results are compared to those obtained using a Polynomial Chaos (PC) surrogate. The possible differences in applicability are highlighted as well. For more on the regression method used here, see PolynomialRegressionTrainer while details of Polynomial Chaos are available under PolynomialChaos.

Problem Statement

The full-order model in this example is essentially the same as the one described in Training a Surrogate Model. It is a one-dimensional heat conduction model:

−kd2Tdx2=q ,x∈[0,L]-k\frac{d^2T}{dx^2} = q \,, \quad x\in[0, L]dTdx∣x=0=0T(x=L)=T∞\begin{aligned} \left.\frac{dT}{dx}\right|_{x=0} &= 0 \\ T(x=L) &= T_{\infty} \end{aligned}

where TT is the temperature, kk is the thermal conductivity, LL is the length of the domain, qq is a heat source and T∞T_\infty is the value of the Dirichlet boundary condition. To make the comparison between different surrogate models easier, only the maximum temperature is selected to be the Quantity of Interest (QoI):

Tmax⁡=max⁡x∈[0,L]T(x).\begin{aligned}\\ T_{\max} &= \max_{x\in[0,L]}T(x). \end{aligned}(1)

The problem is parametric in a sense that the solution depends on four input parameters: T=T(x,k,q,L,T∞)T=T(x,k,q,L,T_\infty). Two problem settings are considered in this example. In the first scenario, all of the parameters are assumed to have Uniform distributions (U(a,b)\mathcal{U}(a,b)), while the second considers parameters with Normal distributions (N(μ,σ)\mathcal{N}(\mu,\sigma)). To be more specific the distributions for the two cases are:

ParameterSymbolUniformNormal
Conductivitykk∼U(1,10)\sim\mathcal{U}(1, 10)∼N(5,2)\sim\mathcal{N}(5, 2)
Volumetric Heat Sourceqq∼U(9000,11000)\sim\mathcal{U}(9000, 11000)∼N(10000,500)\sim\mathcal{N}(10000, 500)
Domain SizeLL∼U(0.01,0.05)\sim\mathcal{U}(0.01, 0.05)∼N(0.03,0.01)\sim\mathcal{N}(0.03, 0.01)
Right Boundary TemperatureT∞T_{\infty}∼U(290,310)\sim\mathcal{U}(290, 310)∼N(300,10)\sim\mathcal{N}(300, 10)

The parameters of the uniform distribution are the minimum and maximum bounds, while the parameters of the normal distribution are the mean and standard deviation. It must be mentioned that the maximum temperature can be determined analytically and turns out to be:

Tmax⁡(k,q,L,T∞)=qL22k+T∞.\begin{aligned} T_{\max}(k,q,L,T_{\infty}) &= \frac{qL^2}{2k} + T_{\infty}. \end{aligned}

Using this expression and the previously described probability density functions, the mean (μ\mu) and standard deviation (σ\sigma) of the QoI can be computed for reference:

Table 1: The reference results for the mean and standard deviation of the maximum temperature.

MomentUniformNormal
μTmax\mu_{T_{max}}301.3219301.2547
σTmax\sigma_{T_{max}}5.958510.0011

Solving the problem without uncertain parameters

The first step towards creating a surrogate model is the generation of a full-order model which can solve Eq. (1) with fixed parameter combinations. The complete input file for this case is presented in Listing 1.

Listing 1: Complete input file for the heat equation problem in this study.

[Mesh]
  type = GeneratedMesh
  dim = 1
  nx = 100
  xmax = 1
  elem_type = EDGE3
[]

[Variables]
  [T]
    order = SECOND
    family = LAGRANGE
  []
[]

[Kernels]
  [diffusion]
    type = MatDiffusion
    variable = T
    diffusivity = k
  []
  [source]
    type = BodyForce
    variable = T
    value = 1.0
  []
[]

[Materials]
  [conductivity]
    type = GenericConstantMaterial
    prop_names = k
    prop_values = 2.0
  []
[]

[BCs]
  [right]
    type = DirichletBC
    variable = T
    boundary = right
    value = 300
  []
[]

[Executioner]
  type = Steady
  solve_type = PJFNK
  petsc_options_iname = '-pc_type -pc_hypre_type'
  petsc_options_value = 'hypre boomeramg'
[]

[Postprocessors]
  [max]
    type = NodalExtremeValue
    variable = T
    value_type = max
  []
[]

[Outputs]
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/sub.i)

Training surrogate models

Both surrogate models are constructed using some knowledge about the full-order problem. This means that the full-order problem is solved multiple times with different parameter samples and the value of the QoI is stored from each computation. This step is managed by a master input file which creates parameter samples, transfers them to the sub-application and collects the results from the completed computations. For more information about setting up master input files see Training a Surrogate Model and Parameter Study. The two complete training input files used for the two cases with the two different parameter distributions are available under uniform and normal.

The training phase starts with the definition of the distributions in the Distributions block. The uniform distributions can be defined as:

[Distributions]
  [k_dist]
    type = Uniform
    lower_bound = 1
    upper_bound = 10
  []
  [q_dist]
    type = Uniform
    lower_bound = 9000
    upper_bound = 11000
  []
  [L_dist]
    type = Uniform
    lower_bound = 0.01
    upper_bound = 0.05
  []
  [Tinf_dist]
    type = Uniform
    lower_bound = 290
    upper_bound = 310
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/uniform_train.i)

For the case with normal distributions the block changes to:

[Distributions]
  [k_dist]
    type = Normal
    mean = 5
    standard_deviation = 2
  []
  [q_dist]
    type = Normal
    mean = 10000
    standard_deviation = 500
  []
  [L_dist]
    type = Normal
    mean = 0.03
    standard_deviation = 0.01
  []
  [Tinf_dist]
    type = Normal
    mean = 300
    standard_deviation = 10
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/normal_train.i)

As a next step, several parameter instances are prepared by sampling the underlying distributions. The sampling objects can be defined in the Samplers block. The generation of these parameter samples is different for the two surrogate models. Meanwhile the polynomial chaos uses the samples at specific quadrature points in the parameters space (generated by a QuadratureSampler), the polynomial regression model is trained using samples from a LatinHypercube. It is visible that the number of sample (num_rows) is set in the LatinHypercube to match the number of samples in the tensor-product quadrature set of QuadratureSampler.

[Samplers]
  [pc_sampler]
    type = Quadrature
    order = 8
    distributions = 'k_dist q_dist L_dist Tinf_dist'
    execute_on = PRE_MULTIAPP_SETUP
  []
  [pr_sampler]
    type = LatinHypercube
    distributions = 'k_dist q_dist L_dist Tinf_dist'
    num_rows = 6560
    execute_on = PRE_MULTIAPP_SETUP
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/normal_train.i)

The objects in blocks Controls, MultiApps, Transfers and Reporters are responsible for managing the communication between master and sub-applications, execution of the sub-applications and the collection of the results. For a more detailed description of these blocks see Parameter Study and Training a Surrogate Model.

[MultiApps]
  [pc_sub]
    type = SamplerFullSolveMultiApp
    input_files = sub.i
    sampler = pc_sampler
  []
  [pr_sub]
    type = SamplerFullSolveMultiApp
    input_files = sub.i
    sampler = pr_sampler
  []
[]

[Controls]
  [pc_cmdline]
    type = MultiAppSamplerControl
    multi_app = pc_sub
    sampler = pc_sampler
    param_names = 'Materials/conductivity/prop_values Kernels/source/value Mesh/xmax BCs/right/value'
  []
  [pr_cmdline]
    type = MultiAppSamplerControl
    multi_app = pr_sub
    sampler = pr_sampler
    param_names = 'Materials/conductivity/prop_values Kernels/source/value Mesh/xmax BCs/right/value'
  []
[]

[Transfers]
  [pc_data]
    type = SamplerReporterTransfer
    from_multi_app = pc_sub
    sampler = pc_sampler
    stochastic_reporter = results
    from_reporter = 'max/value'
  []
  [pr_data]
    type = SamplerReporterTransfer
    from_multi_app = pr_sub
    sampler = pr_sampler
    stochastic_reporter = results
    from_reporter = 'max/value'
  []
[]

[Reporters]
  [results]
    type = StochasticReporter
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/normal_train.i)

The next step is to set up two Trainer objects to generate the surrogate models from the available data. This can be done in the Trainers block. It is visible that both examples use the data from Sampler and Reporter objects. A polynomial chaos surrogate of order 8 and a polynomial regression surrogate with a polynomial of degree at most 4 is used in this study. The PolynomialChaosTrainer also needs knowledge about the underlying parameter distributions to be able to select matching polynomials.

[Trainers]
  [pc_max]
    type = PolynomialChaosTrainer
    execute_on = timestep_end
    order = 8
    distributions = 'k_dist q_dist L_dist Tinf_dist'
    sampler = pc_sampler
    response = results/pc_data:max:value
  []
  [pr_max]
    type = PolynomialRegressionTrainer
    regression_type = "ols"
    execute_on = timestep_end
    max_degree = 4
    sampler = pr_sampler
    response = results/pr_data:max:value
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/uniform_train.i)

As a last step in the training process, the important parameters of the trained surrogates are saved into .rd files. These files can be used to construct the surrogate models again without the need to carry out the training process from the beginning.

[Outputs]
  [pc_out]
    type = SurrogateTrainerOutput
    trainers = 'pc_max'
    execute_on = FINAL
  []
  [pr_out]
    type = SurrogateTrainerOutput
    trainers = 'pr_max'
    execute_on = FINAL
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/normal_train.i)

Evaluation of surrogate models

To evaluate surrogate models, a new master input file has to be created for uniform and normal parameter distributions. The input files contain testing distributions for the parameters defined in the Distributions block. In this study, the training distributions are used for the testing of the surrogates as well. Both surrogate models are tested using the same parameter samples. These samples are selected using LatinHypercube defined in the Samplers block. Since the surrogate models are orders of magnitude faster than the full-order model, 100,000100,000 samples are selected for testing (compared to 6,5606,560 used for training).

[Samplers]
  [sample]
    type = LatinHypercube
    num_rows = 100000
    distributions = 'k_dist q_dist L_dist Tinf_dist'
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/normal_surr.i)

As a next step, two object are created in the Surrogates block for the two surrogate modeling techniques. Both of them are constructed using the information available within the corresponding .rd files.

[Surrogates]
  [pc_max]
    type = PolynomialChaos
    filename = 'normal_train_pc_out_pc_max.rd'
  []
  [pr_max]
    type = PolynomialRegressionSurrogate
    filename = 'normal_train_pr_out_pr_max.rd'
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/normal_surr.i)

These surrogate models can be evaluated at the points defined in the testing sample batch. This is done using objects in the Reporters block.

[Reporters]
  [pc_max_res]
    type = EvaluateSurrogate
    model = pc_max
    sampler = sample
    parallel_type = ROOT
  []
  [pr_max_res]
    type = EvaluateSurrogate
    model = pr_max
    sampler = sample
    parallel_type = ROOT
  []
  [pr_max_stats]
    type = StatisticsReporter
    reporters = 'pr_max_res/pr_max'
    compute = 'mean stddev'
  []
  [pc_max_stats]
    type = PolynomialChaosReporter
    pc_name = 'pc_max'
    statistics = 'mean stddev'
  []
[]
(contrib/moose/modules/stochastic_tools/examples/surrogates/polynomial_regression/normal_surr.i)

Results and Analysis

In this section the results from the different surrogate models are provided. They are compared to the reference results summarized in Table 1. A short analysis of the results is provided as well to showcase potential issues the user might encounter when using polynomial regression.

Uniform parameter distributions

First, the case with parameters having uniform distributions are investigated. The statistical moments obtained by the execution of the surrogate model are summarized in Table 2.

Table 2: Comparison of the statistical moments from different surrogate models assuming uniform parameter distributions.

MomentReferencePoly. ChaosPoly. Reg. (deg. 4)Poly. Reg. (deg. 8)
μTmax\mu_{T_{max}}301.3219301.3218301.3234301.3220
σTmax\sigma_{T_{max}}5.95855.95865.96255.9655

It can be observed that the polynomial chaos surrogate gives results closer to the reference values. It is also visible that by increasing the polynomial order for the regression, the accuracy in the standard deviation slightly decreases. The histogram of the results is presented in Figure 1. It is important to mention that the results for the polynomial regression surrogate were obtained using max_degree=4. It is apparent that the two methods give similar solutions.

Figure 1: Histogram of the maximum temperature coming from the Monte Carlo run using the surrogate models and assuming uniform parameter distributions.

Normal parameter distributions

Next, the case with normally distributed parameters is analyzed. The statistical moments of the results from testing the surrogate model are summarized in Table 3.

Table 3: Comparison of the statistical moments from different surrogate models assuming normal distributions.

MomentReferencePoly. ChaosPoly. Reg. (deg. 4)Poly. Reg. (deg. 8)
μTmax\mu_{T_{max}}301.2547301.3162301.5663301.5810
σTmax\sigma_{T_{max}}10.001110.112511.291230.1675

It is visible that polynomial chaos surrogate gives the closest results to the reference values. Furthermore, the increase in the polynomial degree for the regression leads to a decrease in accuracy for both the mean and the standard deviation. This behavior is often referred to as overfitting which decreases the accuracy with the increasing model parameters. The histogram of the results is presented in Figure 2. It is important to mention that the results for the polynomial regression surrogate were obtained using max_degree=4. It is apparent that the two methods give similar solutions, however the tails of the histogram of the polynomial regression are longer.

Figure 2: Histogram of the maximum temperature coming from the Monte Carlo run using the surrogate models and assuming normal parameter distributions.