generatePCEcoeff: Determine the PCE coefficients using a numerical quadrature... In GPC: Generalized Polynomial Chaos

Description

The inner product used to determine the PCE coefficient is approximated using a numerical quadrature where the integral is approximated with a weighted sum of the function evaluations on some quadrature points.

Usage

 `1` ```generatePCEcoeff(m,n,y,PhiZn,weights) ```

Arguments

 `m` Number of PCE terms `n` Number of quadrature points `y` A n-tuple vector containing the solution of the function evaluated at the n quadrature points `PhiZn` A (m x n) matrix containing the m d-dimensional multivariate polynomial evaluated at the n quadrature points. `weights` A n-tuple vector containing the quadrature weights associated with each quadrature point

Value

 `ym` A m-tuple vector containing the PCE coefficients, ordered according to the canonical sequence of the multivariate polynomial expansion

Jordan Ko

References

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

`getM,CreateQuadrature`
 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57``` ```rm(list = setdiff(ls(), lsf.str())) ### CASE 1: ### Ishigami function definition: y= 1 + Phi_1(x1)*Phi_1(x2) + Phi_3(x2) 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 solutions for the Ishigami 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) # random variable definitgeion d <- 3 # number of random variables L <- 4 # quadrature level in each dimention. # could be anisotropic eg c(3,4,5) for full quadrature P <- L-1 # maximum polynomial expansion (cardinal order) M <- getM(d,P) # number of PCE termrs ParamDistrib <- NULL # ParamDistrib<- list(beta=rep(0.0,d),alpha=rep(0.0,d)) # PCE definition QuadType <- "FULL" # type of quadrature QuadPoly <- rep("LEGENDRE",d) # polynomial to use ExpPoly <- rep("LEGENDRE",d) # polynomial to use # QuadType <- "SPARSE" # type of quadrature # QuadPoly <- 'ClenshawCurtis' # polynomial to use # ExpPoly <- rep("LEGENDRE",d) # polynomial to use Quadrature = CreateQuadrature(d,L,QuadPoly,ExpPoly,QuadType,ParamDistrib) # quadrature # function sampling y <- Model(Quadrature\$QuadNodes) # Ishigami function d=3 # generate PCE coefficients PCE = generatePCEcoeff(M,Quadrature\$QuadSize,y,Quadrature\$PolyNodes,Quadrature\$QuadWeights) # PCE # getting mean and variance PCEMean = PCE\$PCEcoeff[1] PCEVar = sum(PCE\$PCEcoeff[2:M]^2*PCE\$PhiIJ[2:M]) # # Sobol' sensitivity analysis Index <- indexCardinal(d,P) Sobol = getSobol(d,Index,PCE\$PCEcoeff,PCE\$PhiIJ) cat('Mean and variance normalized absolute errors are', abs(PCEMean-MeanExact)/MeanExact, 'and', abs(PCEVar-VarExact)/VarExact,'\n') cat('PCE Sobol indices are ',Sobol\$Values,'\n') cat('Sobol absolte errors are ',abs(Sobol\$Values-SobolExact),'\n') ```