R Graphical Manual

Browse All

Last data update: 2014.03.03

R: Solve polynomial chaos expansion with numerical quadrature

GPCE.quad

R Documentation

Solve polynomial chaos expansion with numerical quadrature

Description

The GPCE.quad function implements a polynomial chaos expansion of a given model or an external model. The strategy for the expansion of the model into a polynomial chaos basis is the Gauss quadrature method where the Gauss quadrature points are used to estimate the integrales corresponding to the coefficients of the expansion. A statistical and a global sensitivity analysis of the model is then carried out.

Usage

GPCE.quad(InputDim,PCSpace,InputDistrib,ParamDistrib,
Model,ModelParam,Output,DesignInput,p,ExpPoly,QuadType,QuadPoly,QuadLevel)

Arguments

`InputDim`	Number of input random variables
`PCSpace`	To decide if it should be removed after discussion with Miguel..
`InputDistrib`	Distribution of the input. Options are `Gaussian`, `Uniform`, `Beta`, `Gamma`
`ParamDistrib`	Shape parameters of the random variables.
`Model`	A function defining the model to analyze or NULL if the model is external
`ModelParam`	Optional parameters for function evaluation.
`Output`	The results of the model evaluation at the quadrature points.
`DesignInput`	To decide if it should be removed after discussion with Miguel.
`p`	The order of the polynomial chaos expansion.

order)

`ExpPoly`	The polynomial used in the expansion. Options are `Hermite`, `Legendre`, `Jacobi`, `Laguerre`
`QuadType`	Type of quadrature. Options are `SPARSE` or `FULL`
`QuadPoly`	The type of one-dimensional quadrature rule to be used. For `SPARSE`, one can use `ClenshawCurtis` or `Fejer`. For `FULL`, one could choose `Hermite`,`Legendre`,`Jacobi` or `Laguerre`
`QuadLevel`	Level of quadrature used in the approximation.

Value

`Designs`	A list containing the quadrature size, n, a vector of the n quadrature weights, and a n * d matrix of the quadrature points
`Output`	Vector of the model output
`PCEcoeff`	Matrix of the kept sparse polynomial basis. TruncSet_ij is the jth polynomial degree associated to the ith variable
`Moments`	A list containing the fourth first moments of the output: mean, variance, standard deviation, skweness and kurtosis
`Sensitivity`	A list containing the sobol sensitivity indices (Values) and the normalized Sobol nominal and total sensitivity indices (S and ST)

Author(s)

Jordan Ko

References

J. Ko, D. Lucor and P. Sagaut, 2008, On Sensitivity of Two-Dimensional Spatially Developing Mixing Layer With Respect to Uncertain Inflow Conditions, Physics of Fluids, 20(7), 07710201-07710220.

J. Ko, 2009, Applications of the generalized polynomial chaos to the numerical simulationof stochastic shear flows, Doctoral thesis, University of Paris VI.

J. Ko, D. Lucor and P. Sagaut, 2011, Effects of base flow uncertainty on Couette flow stability, Computers and Fluids, 43(1), 82-89.

Examples


# CASE 1: Validating Legendre PCE with Ishigami function definition
rm(list = setdiff(ls(), lsf.str()))
Model <- function(x){
  a <- 3
  b <- 7
  res <- sin(pi*x[1,]) + a*sin(pi*x[2,])^2 + b*(pi*x[3,])^4*sin(pi*x[1,])
  return(res)
}

# Determine the exact solution for the Ishigami test function 
a <- 3
b <- 7
MeanExact <- a/2
VarExact <- a^2/8 + b*pi^4/5 + b^2*pi^8/18 + 1/2
SobolExact <- c(b*pi^4/5 + b^2*pi^8/50 + 1/2, a^2/8, 0, 0, 8*b^2*pi^8/225, 0, 0)

d = 3
l = 4
ResultObject=GPCE.quad(Model=Model,InputDim=d,PCSpace="Uniform",InputDistrib=rep('Uniform',d),
                       DesignInput=NULL,p=c(l-1),ExpPoly=rep("LEGENDRE",d),QuadType=c("FULL"),
                       QuadPoly=rep("LEGENDRE",d),QuadLevel=c(l),ParamDistrib=NULL,Output=NULL)
names(ResultObject)
cat('Mean and variance normalized absolute errors are', 
  abs(ResultObject$Moments$PCEMean-MeanExact)/MeanExact,
  'and', abs(ResultObject$Moments$PCEVar-VarExact)/VarExact,'\n')
cat('PCE Sobol indices are ',ResultObject$Sensitivity$Values,'\n')
cat('Sobol absolte errors are ',abs(ResultObject$Sensitivity$Values-SobolExact),'\n')

# CASE 2: Validating Hermite PCE with orthogonal Hermite functions
### Model is a R function as a sum of multivariate Hermite polynomials 

Model <- function(x,param){
    d <- param$d
    p <- param$p
    PCETrue <- param$PCETrue
    
    n <- dim(x)[2]
    index <- indexCardinal(d,p)
    PHerm <- hermite.he.polynomials(p, normalized=FALSE)
    y <- rep(0,n)
    
    for (nn in seq(1,n)){
        for (mm in seq(1,getM(d,p))){
            tmp <- 1;
            for (dd in seq(1,d))
            {
            tmp = tmp * unlist(polynomial.values(PHerm[index[dd,mm]+1],x[dd,nn]))
            }     
            y[nn] = y[nn] + PCETrue[mm]*tmp
        }
    }
    return(y)
}

## Problem definition
d = 2;          # random dimension 
l = 3;          # quadrature level
p = l - 1;      # polynomial order of expansion
m = getM(d,p);  # size of polynomial expansion

## Model definition
ModelParam <- NULL
ModelParam$d <- d
ModelParam$p <- p
ModelParam$PCETrue <- sample(seq(1,m),m,replace = FALSE)

## CASE 1: The model is directly evaluated from the GPCE.quad function
ResultObject=GPCE.quad(InputDim=d,PCSpace="Normal",InputDistrib=rep('Gaussian',d),
                       DesignInput=NULL,p=c(p),ExpPoly=rep("HERMITE",d),QuadType=c("FULL"),
                       QuadPoly=rep("HERMITE",d),QuadLevel=c(l),ParamDistrib=NULL,Output=NULL,
                       Model=Model,ModelParam=ModelParam)
cat("The exact PCE coefficients are: \n")
cat(ModelParam$PCETrue,"\n")
cat("The estimated PCE coefficients are: \n")
cat(ResultObject$PCEcoeff,"\n")

## CASE 2: Model is evaluated separately from the GPCE.quad function
# First, the quadrature points are determined from the GPCE.quad function
NoModelResult=GPCE.quad(InputDim=d,PCSpace="Normal",InputDistrib=rep('Gaussian',d),
                       DesignInput=NULL,p=c(p),ExpPoly=rep("HERMITE",d),QuadType=c("FULL"),
                       QuadPoly=rep("HERMITE",d),QuadLevel=c(l),ParamDistrib=NULL,Output=NULL)
cat("The quadrature points can be determined from the Design variable of the output below: \n")
cat(names(NoModelResult),"\n")

# Second, the model is evalauted at the  quadrature points and stored in Output
Output <- Model(NoModelResult$Design$QuadNodes,ModelParam)

# Third, the model output is passed back to GPCE.quad, along with DesignInput and Output
cat("After Design$QuadNodes are evaluated and stored in Output, 
the results is passed back to GPCE.quad:\n")
NoModelResult=GPCE.quad(InputDim=d,PCSpace = "Normal",InputDistrib=rep('Gaussian',d),
                        DesignInput=NoModelResult$Design$QuadNodes,p=c(p),
                        ExpPoly=rep("HERMITE",d),QuadType=c("FULL"),
                        QuadPoly=rep("HERMITE",d),QuadLevel=c(l),
                        ParamDistrib=NULL,Output=Output)
cat("The exact PCE coefficients are:\n")
cat(ModelParam$PCETrue,"\n")
cat("The estimated PCE coefficients are:\n")
cat(NoModelResult$PCEcoeff,"\n")

### Run the algorithm with the lar regression method
# ResultObject=GPCE.lar(Model=Model, PCSpace="Gaussian", InputDim=d, InputDistrib=rep("Gaussian",d))


# CASE 3: Validating Jacobi PCE with orthogonal Jacobi functions 
rm(list = setdiff(ls(), lsf.str()))

# Multivariate Jacobi polynomial test
Model <- function(x){
    a = 2;
    b = 3;
    res <- 2 + 
        3*(2*(a+1)+(a+b+2)*(x[1,]-1))/2 + 
        5*(2*(a+1)+(a+b+2)*(x[2,]-1))/2 + 
        7*(2*(a+1)+(a+b+2)*(x[3,]-1))/2 + 
        11.*(4*(a+1)*(a+2)+4*(a+b+3)*(a+2)*(x[1,]-1)+(a+b+3)*(a+b+4)*(x[1,]-1)^2)/8 + 
        13*(2*(a+1)+(a+b+2)*(x[1,]-1))*(2*(a+1)+(a+b+2)*(x[2,]-1))/4 + 
        17*(2*(a+1)+(a+b+2)*(x[1,]-1))*(2*(a+1)+(a+b+2)*(x[3,]-1))/4 +         
        19.*(4*(a+1)*(a+2)+4*(a+b+3)*(a+2)*(x[2,]-1)+(a+b+3)*(a+b+4)*(x[2,]-1)^2)/8 + 
        23*(2*(a+1)+(a+b+2)*(x[2,]-1))*(2*(a+1)+(a+b+2)*(x[3,]-1))/4 + 
        29.*(4*(a+1)*(a+2)+4*(a+b+3)*(a+2)*(x[3,]-1)+(a+b+3)*(a+b+4)*(x[3,]-1)^2)/8
    return(res)
}

# x is a random variabls following the beta distribution with the pdf
# f(x,a,b) = (1-x)^a(1+x)^b/2^(a+b+1)/Beta(a+1,b+1), where Beta is the beta function 

d = 3
l = 3
ParamDistrib <- NULL
ParamDistrib$alpha <- rep(2,d)  # 1st shape parameters for the random variables
ParamDistrib$beta <- rep(3,d)   # 2nd shape parameters for the random variables

ResultObject=GPCE.quad(Model=Model,InputDim=d,PCSpace="Beta",InputDistrib=rep('Beta',d),
                       DesignInput=NULL,p=c(l-1),ExpPoly=rep("JACOBI",d),QuadType=c("FULL"),
                       QuadPoly=rep("JACOBI",d),QuadLevel=c(l),ParamDistrib,Output=NULL)

cat('Estimated PCE coefficients are', ResultObject$PCEcoeff) 

# CASE 4: Validating Laguerre PCE with orthogonal Laguerre functions 

rm(list = setdiff(ls(), lsf.str()))

# Multivariate Laguerre polynomial test
Model <- function(x){
    a <- 2
    res <- 2 + 
      3*(-x[1,]+a+1) +
      5*(-x[2,]+a+1) + 
      7*(-x[3,]+a+1) +
      11*(x[1,]^2/2-(a+2)*x[1,]+(a+2)*(a+1)/2) + 
      13*(-x[1,]+a+1)*(-x[2,]+a+1) + 
      17*(-x[1,]+a+1)*(-x[3,]+a+1) + 
      19*(x[2,]^2/2-(a+2)*x[2,]+(a+2)*(a+1)/2) + 
      21*(-x[2,]+a+1)*(-x[3,]+a+1) +
      29*(x[3,]^2/2-(a+2)*x[3,]+(a+2)*(a+1)/2) 
    return(res)
}

# x is a random variabls following the gamma distribution with scale parameter theta = 1 and shape 
# parameter k = a + 1 giving the pdf f(x,a) = x^(k-1)*exp(-x)/Gamma(k) = x^a*exp(-x)/Gamma(a+1), 
# where Gamma is the gamma function 

a = 2
d = 3
l = 3

ParamDistrib = NULL
ParamDistrib$alpha = rep(a,d)

ResultObject=GPCE.quad(Model=Model,InputDim=d,PCSpace="Gamma",InputDistrib=rep('Gamma',d),
                       DesignInput=NULL,p=c(l-1),ExpPoly=rep("LAGUERRE",d),QuadType=c("FULL"),
                       QuadPoly=rep("LAGUERRE",d),QuadLevel=c(l),ParamDistrib,Output=NULL)

cat('Estimated PCE coefficients are', ResultObject$PCEcoeff)

Results


R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(GPC)
Loading required package: randtoolbox
Loading required package: rngWELL
This is randtoolbox. For overview, type 'help("randtoolbox")'.
Loading required package: orthopolynom
Loading required package: polynom
Loading required package: ks
Loading required package: KernSmooth
KernSmooth 2.23 loaded
Copyright M. P. Wand 1997-2009
Loading required package: misc3d
Loading required package: mvtnorm
Loading required package: rgl
Loading required package: lars
Loaded lars 1.2

> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/GPC/GPCE.quad.Rd_%03d_medium.png", width=480, height=480)
> ### Name: GPCE.quad
> ### Title: Solve polynomial chaos expansion with numerical quadrature
> ### Aliases: GPCE.quad
> 
> ### ** Examples
> 
> 
> # CASE 1: Validating Legendre PCE with Ishigami function definition
> rm(list = setdiff(ls(), lsf.str()))
> Model <- function(x){
+   a <- 3
+   b <- 7
+   res <- sin(pi*x[1,]) + a*sin(pi*x[2,])^2 + b*(pi*x[3,])^4*sin(pi*x[1,])
+   return(res)
+ }
> 
> # Determine the exact solution for the Ishigami test function 
> a <- 3
> b <- 7
> MeanExact <- a/2
> VarExact <- a^2/8 + b*pi^4/5 + b^2*pi^8/18 + 1/2
> SobolExact <- c(b*pi^4/5 + b^2*pi^8/50 + 1/2, a^2/8, 0, 0, 8*b^2*pi^8/225, 0, 0)
> 
> d = 3
> l = 4
> ResultObject=GPCE.quad(Model=Model,InputDim=d,PCSpace="Uniform",InputDistrib=rep('Uniform',d),
+                        DesignInput=NULL,p=c(l-1),ExpPoly=rep("LEGENDRE",d),QuadType=c("FULL"),
+                        QuadPoly=rep("LEGENDRE",d),QuadLevel=c(l),ParamDistrib=NULL,Output=NULL)
Quadrature level assumed to be isotropic.
> names(ResultObject)
[1] "PCEcoeff"    "PhiIJ"       "Moments"     "Sensitivity" "Args"       
[6] "Output"      "Design"     
> cat('Mean and variance normalized absolute errors are', 
+   abs(ResultObject$Moments$PCEMean-MeanExact)/MeanExact,
+   'and', abs(ResultObject$Moments$PCEVar-VarExact)/VarExact,'\n')
Mean and variance normalized absolute errors are 0.1257323 and 0.2297988 
> cat('PCE Sobol indices are ',ResultObject$Sensitivity$Values,'\n')
PCE Sobol indices are  10622 0.709061 1.710349e-27 8.9161e-28 9377.792 7.511127e-29 1.331203e-30 
> cat('Sobol absolte errors are ',abs(ResultObject$Sensitivity$Values-SobolExact),'\n')
Sobol absolte errors are  1186.364 0.415939 1.710349e-27 8.9161e-28 7153.338 7.511127e-29 1.331203e-30 
> 
> # CASE 2: Validating Hermite PCE with orthogonal Hermite functions
> ### Model is a R function as a sum of multivariate Hermite polynomials 
> 
> Model <- function(x,param){
+     d <- param$d
+     p <- param$p
+     PCETrue <- param$PCETrue
+     
+     n <- dim(x)[2]
+     index <- indexCardinal(d,p)
+     PHerm <- hermite.he.polynomials(p, normalized=FALSE)
+     y <- rep(0,n)
+     
+     for (nn in seq(1,n)){
+         for (mm in seq(1,getM(d,p))){
+             tmp <- 1;
+             for (dd in seq(1,d))
+             {
+             tmp = tmp * unlist(polynomial.values(PHerm[index[dd,mm]+1],x[dd,nn]))
+             }     
+             y[nn] = y[nn] + PCETrue[mm]*tmp
+         }
+     }
+     return(y)
+ }
> 
> ## Problem definition
> d = 2;          # random dimension 
> l = 3;          # quadrature level
> p = l - 1;      # polynomial order of expansion
> m = getM(d,p);  # size of polynomial expansion
> 
> ## Model definition
> ModelParam <- NULL
> ModelParam$d <- d
> ModelParam$p <- p
> ModelParam$PCETrue <- sample(seq(1,m),m,replace = FALSE)
> 
> ## CASE 1: The model is directly evaluated from the GPCE.quad function
> ResultObject=GPCE.quad(InputDim=d,PCSpace="Normal",InputDistrib=rep('Gaussian',d),
+                        DesignInput=NULL,p=c(p),ExpPoly=rep("HERMITE",d),QuadType=c("FULL"),
+                        QuadPoly=rep("HERMITE",d),QuadLevel=c(l),ParamDistrib=NULL,Output=NULL,
+                        Model=Model,ModelParam=ModelParam)
Quadrature level assumed to be isotropic.
> cat("The exact PCE coefficients are: \n")
The exact PCE coefficients are: 
> cat(ModelParam$PCETrue,"\n")
6 4 1 5 2 3 
> cat("The estimated PCE coefficients are: \n")
The estimated PCE coefficients are: 
> cat(ResultObject$PCEcoeff,"\n")
6 4 1 5 2 3 
> 
> ## CASE 2: Model is evaluated separately from the GPCE.quad function
> # First, the quadrature points are determined from the GPCE.quad function
> NoModelResult=GPCE.quad(InputDim=d,PCSpace="Normal",InputDistrib=rep('Gaussian',d),
+                        DesignInput=NULL,p=c(p),ExpPoly=rep("HERMITE",d),QuadType=c("FULL"),
+                        QuadPoly=rep("HERMITE",d),QuadLevel=c(l),ParamDistrib=NULL,Output=NULL)
Quadrature level assumed to be isotropic.
> cat("The quadrature points can be determined from the Design variable of the output below: \n")
The quadrature points can be determined from the Design variable of the output below: 
> cat(names(NoModelResult),"\n")
Design Design2Eval Args 
> 
> # Second, the model is evalauted at the  quadrature points and stored in Output
> Output <- Model(NoModelResult$Design$QuadNodes,ModelParam)
> 
> # Third, the model output is passed back to GPCE.quad, along with DesignInput and Output
> cat("After Design$QuadNodes are evaluated and stored in Output, 
+ the results is passed back to GPCE.quad:\n")
After Design$QuadNodes are evaluated and stored in Output, 
the results is passed back to GPCE.quad:
> NoModelResult=GPCE.quad(InputDim=d,PCSpace = "Normal",InputDistrib=rep('Gaussian',d),
+                         DesignInput=NoModelResult$Design$QuadNodes,p=c(p),
+                         ExpPoly=rep("HERMITE",d),QuadType=c("FULL"),
+                         QuadPoly=rep("HERMITE",d),QuadLevel=c(l),
+                         ParamDistrib=NULL,Output=Output)
Quadrature level assumed to be isotropic.
Warning message:
In GPCE.quad(InputDim = d, PCSpace = "Normal", InputDistrib = rep("Gaussian",  :
  Outputs evalauted at the correct quadrature points. 

> cat("The exact PCE coefficients are:\n")
The exact PCE coefficients are:
> cat(ModelParam$PCETrue,"\n")
6 4 1 5 2 3 
> cat("The estimated PCE coefficients are:\n")
The estimated PCE coefficients are:
> cat(NoModelResult$PCEcoeff,"\n")
6 4 1 5 2 3 
> 
> ### Run the algorithm with the lar regression method
> # ResultObject=GPCE.lar(Model=Model, PCSpace="Gaussian", InputDim=d, InputDistrib=rep("Gaussian",d))
> 
> 
> # CASE 3: Validating Jacobi PCE with orthogonal Jacobi functions 
> rm(list = setdiff(ls(), lsf.str()))
> 
> # Multivariate Jacobi polynomial test
> Model <- function(x){
+     a = 2;
+     b = 3;
+     res <- 2 + 
+         3*(2*(a+1)+(a+b+2)*(x[1,]-1))/2 + 
+         5*(2*(a+1)+(a+b+2)*(x[2,]-1))/2 + 
+         7*(2*(a+1)+(a+b+2)*(x[3,]-1))/2 + 
+         11.*(4*(a+1)*(a+2)+4*(a+b+3)*(a+2)*(x[1,]-1)+(a+b+3)*(a+b+4)*(x[1,]-1)^2)/8 + 
+         13*(2*(a+1)+(a+b+2)*(x[1,]-1))*(2*(a+1)+(a+b+2)*(x[2,]-1))/4 + 
+         17*(2*(a+1)+(a+b+2)*(x[1,]-1))*(2*(a+1)+(a+b+2)*(x[3,]-1))/4 +         
+         19.*(4*(a+1)*(a+2)+4*(a+b+3)*(a+2)*(x[2,]-1)+(a+b+3)*(a+b+4)*(x[2,]-1)^2)/8 + 
+         23*(2*(a+1)+(a+b+2)*(x[2,]-1))*(2*(a+1)+(a+b+2)*(x[3,]-1))/4 + 
+         29.*(4*(a+1)*(a+2)+4*(a+b+3)*(a+2)*(x[3,]-1)+(a+b+3)*(a+b+4)*(x[3,]-1)^2)/8
+     return(res)
+ }
> 
> # x is a random variabls following the beta distribution with the pdf
> # f(x,a,b) = (1-x)^a(1+x)^b/2^(a+b+1)/Beta(a+1,b+1), where Beta is the beta function 
> 
> d = 3
> l = 3
> ParamDistrib <- NULL
> ParamDistrib$alpha <- rep(2,d)  # 1st shape parameters for the random variables
> ParamDistrib$beta <- rep(3,d)   # 2nd shape parameters for the random variables
> 
> ResultObject=GPCE.quad(Model=Model,InputDim=d,PCSpace="Beta",InputDistrib=rep('Beta',d),
+                        DesignInput=NULL,p=c(l-1),ExpPoly=rep("JACOBI",d),QuadType=c("FULL"),
+                        QuadPoly=rep("JACOBI",d),QuadLevel=c(l),ParamDistrib,Output=NULL)
Quadrature level assumed to be isotropic.
> 
> cat('Estimated PCE coefficients are', ResultObject$PCEcoeff) 
Estimated PCE coefficients are 2 3 5 7 11 13 17 19 23 29> 
> # CASE 4: Validating Laguerre PCE with orthogonal Laguerre functions 
> 
> rm(list = setdiff(ls(), lsf.str()))
> 
> # Multivariate Laguerre polynomial test
> Model <- function(x){
+     a <- 2
+     res <- 2 + 
+       3*(-x[1,]+a+1) +
+       5*(-x[2,]+a+1) + 
+       7*(-x[3,]+a+1) +
+       11*(x[1,]^2/2-(a+2)*x[1,]+(a+2)*(a+1)/2) + 
+       13*(-x[1,]+a+1)*(-x[2,]+a+1) + 
+       17*(-x[1,]+a+1)*(-x[3,]+a+1) + 
+       19*(x[2,]^2/2-(a+2)*x[2,]+(a+2)*(a+1)/2) + 
+       21*(-x[2,]+a+1)*(-x[3,]+a+1) +
+       29*(x[3,]^2/2-(a+2)*x[3,]+(a+2)*(a+1)/2) 
+     return(res)
+ }
> 
> # x is a random variabls following the gamma distribution with scale parameter theta = 1 and shape 
> # parameter k = a + 1 giving the pdf f(x,a) = x^(k-1)*exp(-x)/Gamma(k) = x^a*exp(-x)/Gamma(a+1), 
> # where Gamma is the gamma function 
> 
> a = 2
> d = 3
> l = 3
> 
> ParamDistrib = NULL
> ParamDistrib$alpha = rep(a,d)
> 
> ResultObject=GPCE.quad(Model=Model,InputDim=d,PCSpace="Gamma",InputDistrib=rep('Gamma',d),
+                        DesignInput=NULL,p=c(l-1),ExpPoly=rep("LAGUERRE",d),QuadType=c("FULL"),
+                        QuadPoly=rep("LAGUERRE",d),QuadLevel=c(l),ParamDistrib,Output=NULL)
Quadrature level assumed to be isotropic.
> 
> cat('Estimated PCE coefficients are', ResultObject$PCEcoeff) 
Estimated PCE coefficients are 2 3 5 7 11 13 17 19 21 29> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>