PC-KRIGING METAMODELING: MULTIPLE OUTPUTS

This example showcases an application of polynomial chaos-Kriging (PC-Kriging) to the metamodeling of a simply supported beam model with multiple outputs. The model computes the deflections at several points along the length of the beam subjected to a uniform random load.

Package imports

from uqpylab import sessions, display_util
import numpy as np 
import matplotlib.pyplot as plt

Initialize common plotting parameters

display_util.load_plt_defaults()
uq_colors = display_util.get_uq_color_order()

Start a remote UQCloud session

# Start the session
mySession = sessions.cloud()
# (Optional) Get a convenient handle to the command line interface
uq = mySession.cli
# Reset the session
mySession.reset()

Processing .

.

 done!

 uqpylab.sessions :: INFO     :: This is UQ[py]Lab, version 1.00, running on https://uqcloud.ethz.ch. 
                                 UQ[py]Lab is free software, published under the open source BSD 3-clause license.
                                 To request special permissions, please contact:
                                  - Stefano Marelli (marelli@ibk.baug.ethz.ch).
                                 A new session (a5637836a6034cfd88d0a029f5fb4e75) started.

 uqpylab.sessions :: INFO     :: Reset successful.

Set the random seed for reproducibility

uq.rng(100,'twister');

Computational model

The simply supported beam model is shown in the following figure:

The (negative) deflection of the beam at any longitudinal coordinate $s$ is given by:

$$V(s) = -\frac{p \,s (L^3 - 2\, s^2 L + s^3) }{2E b h^3}$$

This computation is carried out by the function simply_supported_beam_9points. The function evaluates the inputs gathered in the $N \times M$ matrix X, where $N$ and $M$ are the numbers of realizations and inputs, respectively. The inputs are given in the following order:

• $b$: beam width $(m)$
• $h$: beam height $(m)$
• $L$: beam length $(m)$
• $E$: Young's modulus $(Pa)$
• $p$: uniform load $(N/m)$

The function returns the beam deflection $V(s_i)$ at nine equally-spaced points along the length $s_i = i \cdot L/10, \; i=1,\ldots,9.$

Create a MODEL object from the simply_supported_beam_9points function:

ModelOpts = {
    'Type': 'Model',
    'ModelFun': 'simply_supported_beam_9points.model',
    'isVectorized': 'true'
}

myModel = uq.createModel(ModelOpts)

Probabilistic input model

The simply supported beam model has five independent input parameters modeled by lognormal random variables. The parameters of the distributions are given in the following table:

Variable	Description	Distribution	Mean	Std. deviation
b	Beam width	Lognormal	0.15 m	7.5 mm
h	Beam height	Lognormal	0.3 m	15 mm
L	Length	Lognormal	5 m	50 mm
E	Young modulus	Lognormal	30000 MPa	4500 MPa
p	Uniform load	Lognormal	10 kN/m	2 kN/m

Define an INPUT object with the following marginals:

InputOpts = {
    'Marginals': [
        {
        'Name': 'b', # beam width
        'Type': 'Lognormal',
        'Moments': [0.15, 0.0075] # (m)
        },
        {
        'Name': 'h', # beam height
        'Type': 'Lognormal',
        'Moments': [0.3, 0.015] # (m)
        },
        {
        'Name': 'L', # beam length
        'Type': 'Lognormal',
        'Moments': [5, 0.05] # (m)
        },
        {
        'Name': 'E', # Young's modulus
        'Type': 'Lognormal',
        'Moments': [3e10, 4.5e9] # (Pa)
        },
        {
        'Name': 'p', # uniform load
        'Type': 'Lognormal',
        'Moments': [1e4, 2e3] # (N/m) # in Kriging example, there is [1e4 1e3]
        }]
}

myInput = uq.createInput(InputOpts)

PC-Kriging metamodel

Select the PCK for metamodeling tool:

MetaOpts = {
    "Type": "Metamodel",
    "MetaType": "PCK"
}

Specify the sampling strategy and the number of sample points for the experimental design:

MetaOpts["ExpDesign"] = {
    "Sampling" : "LHS",
    "NSamples" : 100
}

Set the maximum degree for the polynomial chaos expansion (PCE) trend to $3$:

MetaOpts["PCE"] = {
    "Degree": 3
}

Create the PC-Kriging metamodel:

myPCK = uq.createModel(MetaOpts)

 uqpylab.sessions :: INFO     :: Received intermediate compute request, function: simply_supported_beam_9points.model.

 uqpylab.sessions :: INFO     :: Carrying out local computation...

 uqpylab.sessions :: INFO     :: Local computation complete.

 uqpylab.sessions :: INFO     :: Starting transmission of intermediate compute results ((100, 9))...

 uqpylab.sessions :: INFO     :: Intermediate compute results sent.

Processing .

.

 done!

Validation

The deflections $V(s_i)$ at the nine points are plotted for three realizations of the random inputs. Relative length units are used for comparison, because $L$ is one of the random inputs.

Create and evaluate a validation sample

Generate a validation sample:

Nval = 3
Xval = uq.getSample(myInput, Nval)

Evaluate the original computational model at the validation points:

Yval = uq.evalModel(myModel, Xval)

Evaluate the PC-Kriging metamodel at the same validation set points:

YPCK = uq.evalModel(myPCK, Xval)

Visualize results

For each sample point of the validation set $\mathbf{x}^{(i)}$, the simply supported beam deflection $\mathcal{M}(\mathbf{x}^{(i)})$ is compared against the one predicted by the Kriging metamodel $\mathcal{M}^{PCK}(\mathbf{x}^{(i)})$ (its mean):

li = np.linspace(0,1,11,endpoint=True)
legend_text = []

for idx in np.arange(0,Nval):
    x=li
    y=np.concatenate([[0],YPCK[idx,:],[0]])
    pp = plt.plot(x,y,':',linewidth=1.0, c=uq_colors[idx])
    # legend_text.append('$\mathrm{\mathcal{M}}(\mathbf{x}^{('+str(idx+1)+')}$)')
    legend_text.append('$M(x^{('+str(idx+1)+')}$)')

    y=np.concatenate([[0],Yval[idx,:],[0]])
    plt.plot(x,y,'-',linewidth=1.0,c=uq_colors[idx])
    # Add text for the plot legend
    legend_text.append('$M^{PCK}(x^{('+str(idx+1)+')}$)')

plt.xlabel('$L_{rel}$')
plt.ylabel('$V(m)$')
plt.grid(True)
plt.legend(legend_text);
plt.show()

Terminate the remote UQCloud session

mySession.quit()

 uqpylab.sessions :: INFO     :: Session a5637836a6034cfd88d0a029f5fb4e75 terminated.

True