Stochastic Polynomial approximation of PDEs with random coefficients - - PowerPoint PPT Presentation

stochastic polynomial approximation of pdes with random
SMART_READER_LITE
LIVE PREVIEW

Stochastic Polynomial approximation of PDEs with random coefficients - - PowerPoint PPT Presentation

Stochastic Polynomial approximation of PDEs with random coefficients Fabio Nobile CSQI-MATHICSE, EPFL, Switzerland and MOX, Politecnico di Milano, Italy Joint work with : R. Tempone, E. von Schwerin (KAUST) L. Tamellini, G. Migliorati (MOX,


slide-1
SLIDE 1

Stochastic Polynomial approximation of PDEs with random coefficients

Fabio Nobile

CSQI-MATHICSE, EPFL, Switzerland and MOX, Politecnico di Milano, Italy

Joint work with: R. Tempone, E. von Schwerin (KAUST)

  • L. Tamellini, G. Migliorati (MOX, Politecnico Milano), J. Beck (UCL)

WS: “Numer. Anal. of Multiscale Problems & Stochastic Modelling” RICAM, Linz, December 12-16, 2011

Italian project FIRB-IDEAS (’09) Advanced Numerical Tech- niques for Uncertainty Quantification in Engineering and Life Science Problems

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 1

slide-2
SLIDE 2

Outline

1

Elliptic PDE with random coefficients

2

Stochastic multivariate polynomial approximation Galerkin projection Collocation on sparse grids Optimized algorithms Numerical results

3

Polynomial approximation by discrete projection on random points

4

Conclusions

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 3

slide-3
SLIDE 3

Elliptic PDE with random coefficients

Outline

1

Elliptic PDE with random coefficients

2

Stochastic multivariate polynomial approximation Galerkin projection Collocation on sparse grids Optimized algorithms Numerical results

3

Polynomial approximation by discrete projection on random points

4

Conclusions

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 4

slide-4
SLIDE 4

Elliptic PDE with random coefficients

Elliptic PDE with random coefficients

Let (Ω, F, P) be a complete probability space. Consider L(u) = F ⇔

  • − div(a(ω, x)∇u(ω, x)) = f (x)

x ∈ D, ω ∈ Ω, u(ω, x) = 0 x ∈ ∂D, ω ∈ Ω where a : Ω × D → R is a second order random field such that a ≥ amin > 0 a.e. in Ω × D f ∈ L2(D) deterministic (could be stochastic as well,

f ∈ L2

P(Ω) ⊗ L2(D))

By Lax-Milgram lemma, for a.e. ω ∈ Ω, there exists a unique solution u(ω, ·) ∈ V ≡ H1

0(D),

and uV ⊗L2

P ≤ CP

amin f L2(D) The uniform coercivity assumption can be relaxed to

E[amin(ω)−p] < ∞ for all p > 0. Useful for lognormal random fields

(see e.g. [Babuska-N.-Tempone ’07, Garvis-Sarkis ’09, Gittelson ’10, Charrier ’11]).

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 5

slide-5
SLIDE 5

Elliptic PDE with random coefficients

Elliptic PDE with random coefficients

Let (Ω, F, P) be a complete probability space. Consider L(u) = F ⇔

  • − div(a(ω, x)∇u(ω, x)) = f (x)

x ∈ D, ω ∈ Ω, u(ω, x) = 0 x ∈ ∂D, ω ∈ Ω where a : Ω × D → R is a second order random field such that a ≥ amin > 0 a.e. in Ω × D f ∈ L2(D) deterministic (could be stochastic as well,

f ∈ L2

P(Ω) ⊗ L2(D))

By Lax-Milgram lemma, for a.e. ω ∈ Ω, there exists a unique solution u(ω, ·) ∈ V ≡ H1

0(D),

and uV ⊗L2

P ≤ CP

amin f L2(D) The uniform coercivity assumption can be relaxed to

E[amin(ω)−p] < ∞ for all p > 0. Useful for lognormal random fields

(see e.g. [Babuska-N.-Tempone ’07, Garvis-Sarkis ’09, Gittelson ’10, Charrier ’11]).

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 5

slide-6
SLIDE 6

Elliptic PDE with random coefficients

Elliptic PDE with random coefficients

Let (Ω, F, P) be a complete probability space. Consider L(u) = F ⇔

  • − div(a(ω, x)∇u(ω, x)) = f (x)

x ∈ D, ω ∈ Ω, u(ω, x) = 0 x ∈ ∂D, ω ∈ Ω where a : Ω × D → R is a second order random field such that a ≥ amin > 0 a.e. in Ω × D f ∈ L2(D) deterministic (could be stochastic as well,

f ∈ L2

P(Ω) ⊗ L2(D))

By Lax-Milgram lemma, for a.e. ω ∈ Ω, there exists a unique solution u(ω, ·) ∈ V ≡ H1

0(D),

and uV ⊗L2

P ≤ CP

amin f L2(D) The uniform coercivity assumption can be relaxed to

E[amin(ω)−p] < ∞ for all p > 0. Useful for lognormal random fields

(see e.g. [Babuska-N.-Tempone ’07, Garvis-Sarkis ’09, Gittelson ’10, Charrier ’11]).

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 5

slide-7
SLIDE 7

Elliptic PDE with random coefficients

Assumption – random field parametrization

The stochastic coefficient a(ω, x) can be parametrized by a finite number of independent random variables a(ω, x) = a(y1(ω), . . . , yN(ω), x) We assume that y has a joint probability density function ρ(y) = N

n=1 ρn(yn) : ΓN → R+

Then u(ω, x) = u(y1(ω), . . . , yN(ω), x) is a deterministic function of the random vector y. Extensions to the case of infinitely many (countable) random variables is possible provided the solution u(y1, . . . , yn, . . . , x) has suitable decay properties w.r.t. n.

(see e.g. [Cohen-Devore-Schwab, ’09,’10])

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 6

slide-8
SLIDE 8

Elliptic PDE with random coefficients

Assumption – random field parametrization

The stochastic coefficient a(ω, x) can be parametrized by a finite number of independent random variables a(ω, x) = a(y1(ω), . . . , yN(ω), x) We assume that y has a joint probability density function ρ(y) = N

n=1 ρn(yn) : ΓN → R+

Then u(ω, x) = u(y1(ω), . . . , yN(ω), x) is a deterministic function of the random vector y. Extensions to the case of infinitely many (countable) random variables is possible provided the solution u(y1, . . . , yn, . . . , x) has suitable decay properties w.r.t. n.

(see e.g. [Cohen-Devore-Schwab, ’09,’10])

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 6

slide-9
SLIDE 9

Elliptic PDE with random coefficients

Assumption – random field parametrization

The stochastic coefficient a(ω, x) can be parametrized by a finite number of independent random variables a(ω, x) = a(y1(ω), . . . , yN(ω), x) We assume that y has a joint probability density function ρ(y) = N

n=1 ρn(yn) : ΓN → R+

Then u(ω, x) = u(y1(ω), . . . , yN(ω), x) is a deterministic function of the random vector y. Extensions to the case of infinitely many (countable) random variables is possible provided the solution u(y1, . . . , yn, . . . , x) has suitable decay properties w.r.t. n.

(see e.g. [Cohen-Devore-Schwab, ’09,’10])

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 6

slide-10
SLIDE 10

Elliptic PDE with random coefficients

Examples

Material with inclusions of random conductivity

a(ω, x) = a0 + N

n=N yn(ω)✶Ωn(x)

with yn ∼ uniform, lognormal, ...

Random, spatially correlated, material properties

a(ω, x) is ∞-dimensional random field (e.g. lognormal), suitably truncated a(ω, x) ≈ amin + e

N

n=1 yn(ω)bn(x)

with yn ∼ N(0, 1), i.i.d.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 7

slide-11
SLIDE 11

Elliptic PDE with random coefficients

Examples

Material with inclusions of random conductivity

a(ω, x) = a0 + N

n=N yn(ω)✶Ωn(x)

with yn ∼ uniform, lognormal, ...

Random, spatially correlated, material properties

a(ω, x) is ∞-dimensional random field (e.g. lognormal), suitably truncated a(ω, x) ≈ amin + e

N

n=1 yn(ω)bn(x)

with yn ∼ N(0, 1), i.i.d.

random field with Lc=1/4

0.2 0.4 0.6 0.8 1 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −1 −0.5 0.5 1 1.5 2 2.5

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 7

slide-12
SLIDE 12

Elliptic PDE with random coefficients

Analytic regularity

Theorem [Back-N.-Tamellini-Tempone ’11, Babuska-N.-Tempone ’05, Cohen-DeVore-Schwab ’09/’10] Assume y bounded (e.g. y ∈ ΓN ≡ [−1, 1]N) Let i = (i1, . . . , iN) ∈ NN and r = (r1, . . . , rN) > 0. Set ri =

n r in n .

assume 1

a ∂ia ∂yiL∞(D) ≤ ri uniformly in y

Then ∂iu

∂yi V ≤ C|i|!( 1 log 2r)i

uniformly in y u : ΓN → V is analytic and can be extended analyticially to Σ =

  • z ∈ CN :

N

  • n=1

rn|zn − ˜ yn| < log 2 for some ˜ y ∈ ΓN

  • Better estimates on analyticity region can be obtained by complex analysis.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 8

slide-13
SLIDE 13

Elliptic PDE with random coefficients

Analytic regularity

Theorem [Back-N.-Tamellini-Tempone ’11, Babuska-N.-Tempone ’05, Cohen-DeVore-Schwab ’09/’10] Assume y bounded (e.g. y ∈ ΓN ≡ [−1, 1]N) Let i = (i1, . . . , iN) ∈ NN and r = (r1, . . . , rN) > 0. Set ri =

n r in n .

assume 1

a ∂ia ∂yiL∞(D) ≤ ri uniformly in y

Then ∂iu

∂yi V ≤ C|i|!( 1 log 2r)i

uniformly in y u : ΓN → V is analytic and can be extended analyticially to Σ =

  • z ∈ CN :

N

  • n=1

rn|zn − ˜ yn| < log 2 for some ˜ y ∈ ΓN

  • Better estimates on analyticity region can be obtained by complex analysis.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 8

slide-14
SLIDE 14

Stochastic multivariate polynomial approximation

Outline

1

Elliptic PDE with random coefficients

2

Stochastic multivariate polynomial approximation Galerkin projection Collocation on sparse grids Optimized algorithms Numerical results

3

Polynomial approximation by discrete projection on random points

4

Conclusions

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 10

slide-15
SLIDE 15

Stochastic multivariate polynomial approximation

Stochastic multivariate polynomial approximation

Idea: approximate the response function u(y, ·) by global multivariate

  • polynomials. Since u(y, ·) is analytic, we expect fast convergence.

Let Λ ⊂ NN be an index set of cardinality |Λ| = M, and consider the multivariate polynomial space PΛ(ΓN) = span N

n=1 y pn n ,

with p = (p1, . . . , pN) ∈ Λ

  • Polynomial approximation

find M particular solutions up ∈ V , ∀p ∈ Λ and build uΛ(y, x) =

  • p∈Λ

up(x)y p1

1 y p2 2 · · · y pN N

Compute statistics of functionals J(u), J ∈ H−1(D) of the solution as E[J(u)] ≈ E[J(uΛ)]

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 11

slide-16
SLIDE 16

Stochastic multivariate polynomial approximation

Stochastic multivariate polynomial approximation

Idea: approximate the response function u(y, ·) by global multivariate

  • polynomials. Since u(y, ·) is analytic, we expect fast convergence.

Let Λ ⊂ NN be an index set of cardinality |Λ| = M, and consider the multivariate polynomial space PΛ(ΓN) = span N

n=1 y pn n ,

with p = (p1, . . . , pN) ∈ Λ

  • Polynomial approximation

find M particular solutions up ∈ V , ∀p ∈ Λ and build uΛ(y, x) =

  • p∈Λ

up(x)y p1

1 y p2 2 · · · y pN N

Compute statistics of functionals J(u), J ∈ H−1(D) of the solution as E[J(u)] ≈ E[J(uΛ)]

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 11

slide-17
SLIDE 17

Stochastic multivariate polynomial approximation

Stochastic multivariate polynomial approximation

Idea: approximate the response function u(y, ·) by global multivariate

  • polynomials. Since u(y, ·) is analytic, we expect fast convergence.

Let Λ ⊂ NN be an index set of cardinality |Λ| = M, and consider the multivariate polynomial space PΛ(ΓN) = span N

n=1 y pn n ,

with p = (p1, . . . , pN) ∈ Λ

  • Polynomial approximation

find M particular solutions up ∈ V , ∀p ∈ Λ and build uΛ(y, x) =

  • p∈Λ

up(x)y p1

1 y p2 2 · · · y pN N

Compute statistics of functionals J(u), J ∈ H−1(D) of the solution as E[J(u)] ≈ E[J(uΛ)]

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 11

slide-18
SLIDE 18

Stochastic multivariate polynomial approximation

Examples of pol. spaces: N = 2, p = 16

Tensor product:

pn ≤ w

2 4 6 8 10 12 14 16 2 4 6 8 10 12 14 16

p1 p2

Tensor Product (TP)

2 4 6 8 10 12 14 16 2 4 6 8 10 12 14 16

p1 p2

Total Degree (TD)

Total degree:

  • n pn ≤ w

Hyperbolic cross:

  • n (pn + 1)

≤ w + 1

2 4 6 8 10 12 14 16 2 4 6 8 10 12 14 16

p1 p2

Hyperbolic Cross (HC)

2 4 6 8 10 12 14 16 2 4 6 8 10 12 14 16

p1 p2

Smolyak (SM)

Smolyak:

  • n f (pn) ≤ f (w)

f (p) =      0, p = 0 1, p = 1 ⌈log2(p)⌉, p ≥ 2 Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 12

slide-19
SLIDE 19

Stochastic multivariate polynomial approximation Galerkin projection

Galerkin projection

[Ghanem-Spanos, Karniadakis et al, Matthies-Keese, Schwab-Todor et al., Knio-Le Maˆ ıtre et al,Babuska et al.,. . . ]

Project the equation L(y)(u) = F onto the subspace V ⊗ PΛ(Γ) Let {ψj}M

j=1 be an orthonormal basis w.r.t. the probability density

ρ(y). Expand uΛ(y) on the basis: uΛ(y) = M

j=1 ujψj(y)

Galerkin formulation Find uj ∈ V , j = 1, . . . M s.t. E

  • L(y)(

M

  • j=1

ujψj)ψi

  • = E[Fψi],

i = 1, . . . , M This approach leads to solving M coupled deterministic problems; difficult to assemble and need good preconditioners.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 13

slide-20
SLIDE 20

Stochastic multivariate polynomial approximation Collocation on sparse grids

Collocation approach

[Smolyak ’63, Griebel et al ’98-’03-’04, Barthelmann-Novak-Ritter ’00, Hesthaven-Xiu ’05, N.-Tempone-Webster ’08, Zabaras et al ’07]

1

Choose a set of points y(j) ∈ Γ, j = 1, . . . , M

2

Compute the solutions uj ∈ V : L(y(j))(uj) = F

3

Interpolate the obtained values: uΛ(y) =

M j=1 ujφj(y).

φj ∈ PΛ(Γ): suitable combinations of Lagrange polynomials Always leads to solving M uncoupled deterministic problems The number M of points needed is larger than the dimension M

  • f the polynomial space (Except for tensor product spaces).

Tensor grids are impractical in high dimension (curse of dimensionality)

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 15

slide-21
SLIDE 21

Stochastic multivariate polynomial approximation Collocation on sparse grids

Collocation approach

[Smolyak ’63, Griebel et al ’98-’03-’04, Barthelmann-Novak-Ritter ’00, Hesthaven-Xiu ’05, N.-Tempone-Webster ’08, Zabaras et al ’07]

1

Choose a set of points y(j) ∈ Γ, j = 1, . . . , M

2

Compute the solutions uj ∈ V : L(y(j))(uj) = F

3

Interpolate the obtained values: uΛ(y) =

M j=1 ujφj(y).

φj ∈ PΛ(Γ): suitable combinations of Lagrange polynomials Always leads to solving M uncoupled deterministic problems The number M of points needed is larger than the dimension M

  • f the polynomial space (Except for tensor product spaces).

Tensor grids are impractical in high dimension (curse of dimensionality)

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 15

slide-22
SLIDE 22

Stochastic multivariate polynomial approximation Collocation on sparse grids

Collocation on a (generalized) Sparse Grid

Let i = [i1, . . . , iN] ∈ NN

+ and m(i) : N+ → N+ an increasing function

1

1D polynomial interpolant operators: U m(in)

n

  • n m(in) abscissas.

We use either

Clenshaw-Curtis (extrema on Chebyshev polynomials) Gauss points w.r.t. the weight ρn, assuming that the probability density factorizes as ρ(y) = N

n=1 ρn(yn)

2

Detail operator: ∆m(in)

n

= U m(in)

n

− U m(in−1)

n

, U m(0)

n

= 0.

3

Sparse grid approximation: on an index set Λ ⊂ NN uΛ =

  • i∈Λ

N

  • n=1

∆m(in)

n

[u]

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 16

slide-23
SLIDE 23

Stochastic multivariate polynomial approximation Collocation on sparse grids

By choosing properly the function m and the set Λ one can obtain a polynomial approximation in any given multivariate polynomial space ([Back-N.-Tamellini-Tempone, 2010]) Examples of sparse grids: N = 2,

  • max. polynomial degree p = 16

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 17

slide-24
SLIDE 24

Stochastic multivariate polynomial approximation Collocation on sparse grids

By choosing properly the function m and the set Λ one can obtain a polynomial approximation in any given multivariate polynomial space ([Back-N.-Tamellini-Tempone, 2010]) Examples of sparse grids: N = 2,

  • max. polynomial degree p = 16

−1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1 −1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1

y1 y2

Tensor Product (TP)

−1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1 −1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1

y1 y2

Total Degree (TD)

−1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1 −1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1

y1 y2

Hyperbolic Cross (HC)

−1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1 −1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1

y1 y2

Smolyak Gauss (SM)

−1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1 −1 −0.8 −0.6 −0.4 −0.2 0.2 0.4 0.6 0.8 1

y1 y2

Smolyak CC (SM)

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 17

slide-25
SLIDE 25

Stochastic multivariate polynomial approximation Collocation on sparse grids

A simple numerical test

1D problem:

  • −(a(x, ω)u(x, ω)′)′ = sin(πx),

x ∈ (0, 1), u(0, ω) = u(1, ω) = 0 Diffusion coefficient a(x, ω) = eγ(x,ω): lognormal with either exponential or Gaussian covariance γ(x, ω): stationary Gaussian random field, mean=0; std=1; correlation length: 0.3 Karhunen-Lo` eve expansion

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −3 −2 −1 1 2 3 x log − conductivity

Exponential cov.

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −3 −2 −1 1 2 3 x log − conductivity

Gaussian cov.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 18

slide-26
SLIDE 26

Stochastic multivariate polynomial approximation Collocation on sparse grids

A simple numerical test

Stochastic Collocation approximation: Gauss-Hermite nodes isotropic TD grid (m(i) = i; Λ ≡ {i : |i|1 ≤ w}) measured error E[u] − E[uΛ]L2ρ comparison with Monte Carlo Exponential cov. Gaussian cov.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 19

slide-27
SLIDE 27

Stochastic multivariate polynomial approximation Optimized algorithms

Optimization of polynomial spaces

Galerkin uSG

Λ

=

  • p∈Λ

up(x)ψp(y) find uSG

Λ

by Galerkin projection

  • f the equation on

PΛ = span{ψp, p ∈ Λ}. Collocation uSC

Λ

=

  • i∈Λ
  • n=1,...,N

∆m(in)

n

[u]. Compute uSC

Λ

by collocation on the corresponding sparse grid Question: What is the best index set Λ in both cases?

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 20

slide-28
SLIDE 28

Stochastic multivariate polynomial approximation Optimized algorithms

Optimization of polynomial spaces

Galerkin uSG

Λ

=

  • p∈Λ

up(x)ψp(y) find uSG

Λ

by Galerkin projection

  • f the equation on

PΛ = span{ψp, p ∈ Λ}. Collocation uSC

Λ

=

  • i∈Λ
  • n=1,...,N

∆m(in)

n

[u]. Compute uSC

Λ

by collocation on the corresponding sparse grid Question: What is the best index set Λ in both cases?

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 20

slide-29
SLIDE 29

Stochastic multivariate polynomial approximation Optimized algorithms

Galerkin projection – best M term approximation

Galerkin optimality: u − uSG

Λ V ⊗L2

ρ(Γ) ≤ C

inf

vΛ∈V ⊗PΛ u − vΛV ⊗L2

ρ(Γ)

Let {ψp, p ∈ NN} be the orthonormal basis of multivariate polynomials w.r.t. the denisty ρ(y) = N

n=1 ρn(yn) and vΛ the

truncated expansion of u vΛ =

  • p∈Λ

upψp, up = E[uψp] Parseval’s identity: u − vΛ2

V ⊗L2

ρ(Γ) =

p/ ∈Λ up2 V

Best M terms approximation The optimal index set Λ of cardinality M is the one that contains the M largest Legendre coefficients upV

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 21

slide-30
SLIDE 30

Stochastic multivariate polynomial approximation Optimized algorithms

Galerkin projection – best M term approximation

Galerkin optimality: u − uSG

Λ V ⊗L2

ρ(Γ) ≤ C

inf

vΛ∈V ⊗PΛ u − vΛV ⊗L2

ρ(Γ)

Let {ψp, p ∈ NN} be the orthonormal basis of multivariate polynomials w.r.t. the denisty ρ(y) = N

n=1 ρn(yn) and vΛ the

truncated expansion of u vΛ =

  • p∈Λ

upψp, up = E[uψp] Parseval’s identity: u − vΛ2

V ⊗L2

ρ(Γ) =

p/ ∈Λ up2 V

Best M terms approximation The optimal index set Λ of cardinality M is the one that contains the M largest Legendre coefficients upV

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 21

slide-31
SLIDE 31

Stochastic multivariate polynomial approximation Optimized algorithms

Galerkin projection – best M term approximation

Galerkin optimality: u − uSG

Λ V ⊗L2

ρ(Γ) ≤ C

inf

vΛ∈V ⊗PΛ u − vΛV ⊗L2

ρ(Γ)

Let {ψp, p ∈ NN} be the orthonormal basis of multivariate polynomials w.r.t. the denisty ρ(y) = N

n=1 ρn(yn) and vΛ the

truncated expansion of u vΛ =

  • p∈Λ

upψp, up = E[uψp] Parseval’s identity: u − vΛ2

V ⊗L2

ρ(Γ) =

p/ ∈Λ up2 V

Best M terms approximation The optimal index set Λ of cardinality M is the one that contains the M largest Legendre coefficients upV

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 21

slide-32
SLIDE 32

Stochastic multivariate polynomial approximation Optimized algorithms

Galerkin projection – best M term approximation

Galerkin optimality: u − uSG

Λ V ⊗L2

ρ(Γ) ≤ C

inf

vΛ∈V ⊗PΛ u − vΛV ⊗L2

ρ(Γ)

Let {ψp, p ∈ NN} be the orthonormal basis of multivariate polynomials w.r.t. the denisty ρ(y) = N

n=1 ρn(yn) and vΛ the

truncated expansion of u vΛ =

  • p∈Λ

upψp, up = E[uψp] Parseval’s identity: u − vΛ2

V ⊗L2

ρ(Γ) =

p/ ∈Λ up2 V

Best M terms approximation The optimal index set Λ of cardinality M is the one that contains the M largest Legendre coefficients upV

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 21

slide-33
SLIDE 33

Stochastic multivariate polynomial approximation Optimized algorithms

Estimate of Fourier coefficients

For the diffusion problem with uniform random variables, the solution u(y) is analytic in Γ = [−1, 1]N and the following estimate of the Legendre coefficients holds [Cohen-DeVore-Schwab ’10, Back-N.-Tamellini-Tempone ’11] upV ≤ C0e−

n gnpn |p|!

p! (1) for some gn > 0, with |p| =

n pn, p! = n pn!.

Then the optimal index set of level w is (TD-FC) Λ(w) =

  • p ∈ NN :
  • n

gnpn − log |p|! p! ≤ w

  • In practice, we estimate numerically the rates gn by inexpensive 1D

analyses.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 24

slide-34
SLIDE 34

Stochastic multivariate polynomial approximation Optimized algorithms

Estimate of Fourier coefficients

For the diffusion problem with uniform random variables, the solution u(y) is analytic in Γ = [−1, 1]N and the following estimate of the Legendre coefficients holds [Cohen-DeVore-Schwab ’10, Back-N.-Tamellini-Tempone ’11] upV ≤ C0e−

n gnpn |p|!

p! (1) for some gn > 0, with |p| =

n pn, p! = n pn!.

Then the optimal index set of level w is (TD-FC) Λ(w) =

  • p ∈ NN :
  • n

gnpn − log |p|! p! ≤ w

  • In practice, we estimate numerically the rates gn by inexpensive 1D

analyses.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 24

slide-35
SLIDE 35

Stochastic multivariate polynomial approximation Optimized algorithms

Estimate of Fourier coefficients

For the diffusion problem with uniform random variables, the solution u(y) is analytic in Γ = [−1, 1]N and the following estimate of the Legendre coefficients holds [Cohen-DeVore-Schwab ’10, Back-N.-Tamellini-Tempone ’11] upV ≤ C0e−

n gnpn |p|!

p! (1) for some gn > 0, with |p| =

n pn, p! = n pn!.

Then the optimal index set of level w is (TD-FC) Λ(w) =

  • p ∈ NN :
  • n

gnpn − log |p|! p! ≤ w

  • In practice, we estimate numerically the rates gn by inexpensive 1D

analyses.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 24

slide-36
SLIDE 36

Stochastic multivariate polynomial approximation Optimized algorithms

Numerical tests

We consider the 1D problem

  • −(a(x, y)u(x, y)′)′ = 1

x ∈ D = (0, 1), y ∈ Γ u(0, y) = u(1, y) = 0, y ∈ Γ with several choices of a(x, y) and compute Θ(u) = u( 1

2).

We compare: (Aniso) TD space: Λ(w) =

  • p ∈ NN :
  • n

gnpn ≤ w

  • .

(Aniso) TD-FC space: Λ(w) =

  • p ∈ NN :

N

  • n=1

gnpn − log |p|! p! ≤ w

  • .

The rates gn have been estimated numerically by inexpensive 1D

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 26

slide-37
SLIDE 37

Stochastic multivariate polynomial approximation Optimized algorithms

Numerical tests

We consider the 1D problem

  • −(a(x, y)u(x, y)′)′ = 1

x ∈ D = (0, 1), y ∈ Γ u(0, y) = u(1, y) = 0, y ∈ Γ with several choices of a(x, y) and compute Θ(u) = u( 1

2).

We compare: (Aniso) TD space: Λ(w) =

  • p ∈ NN :
  • n

gnpn ≤ w

  • .

(Aniso) TD-FC space: Λ(w) =

  • p ∈ NN :

N

  • n=1

gnpn − log |p|! p! ≤ w

  • .

The rates gn have been estimated numerically by inexpensive 1D

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 26

slide-38
SLIDE 38

Stochastic multivariate polynomial approximation Optimized algorithms

A numerical check: a(x, y) = 1 + 0.1xy1 + 0.5x2y2

20 40 60 80 100 10

−20

10

−18

10

−16

10

−14

10

−12

10

−10

10

−8

10

−6

10

−4

10

−2

Legendre coeff TD appr TD−FC appr

Figure: Legendre coeffs of Θ(u) in lexicographic order, with TD and TD-FC estimates

20 40 60 80 10

−9

10

−8

10

−7

10

−6

10

−5

10

−4

10

−3

10

−2

10

−1

best M term TD iso−TD TD−FC

Figure: Convergence plot for Θ(u) − Θ(uM)2

L2

ρ(Γ) w.r.t.

M = |Λ| The Legendre coefficients have been computed with a sufficiently high level sparse grids.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 27

slide-39
SLIDE 39

Stochastic multivariate polynomial approximation Optimized algorithms

Optimization of sparse grids

uM = Sm

Λ [u] =

  • i∈Λ

N

  • n=1

∆m(in)

n

[u]. We use a knapsack problem-approach [Griebel-Knapek ’09, Gerstner-Griebel ’03,

Bungartz-Griebel ’04]: for each multiindex i estimate

∆E(i): how much error decreases if i is added to Λ (error contribution) ∆W (i): how much work, i.e. number of evaluations, increases if i is added to Λ (work contribution) Then estimate the profit of each i as P(i) = ∆E(i) ∆W (i) and build the sparse grid using the set Λ of the M indices with the largest profit.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 28

slide-40
SLIDE 40

Stochastic multivariate polynomial approximation Optimized algorithms

Optimization of sparse grids

uM = Sm

Λ [u] =

  • i∈Λ

N

  • n=1

∆m(in)

n

[u]. We use a knapsack problem-approach [Griebel-Knapek ’09, Gerstner-Griebel ’03,

Bungartz-Griebel ’04]: for each multiindex i estimate

∆E(i): how much error decreases if i is added to Λ (error contribution) ∆W (i): how much work, i.e. number of evaluations, increases if i is added to Λ (work contribution) Then estimate the profit of each i as P(i) = ∆E(i) ∆W (i) and build the sparse grid using the set Λ of the M indices with the largest profit.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 28

slide-41
SLIDE 41

Stochastic multivariate polynomial approximation Optimized algorithms

Optimization of sparse grids

uM = Sm

Λ [u] =

  • i∈Λ

N

  • n=1

∆m(in)

n

[u]. We use a knapsack problem-approach [Griebel-Knapek ’09, Gerstner-Griebel ’03,

Bungartz-Griebel ’04]: for each multiindex i estimate

∆E(i): how much error decreases if i is added to Λ (error contribution) ∆W (i): how much work, i.e. number of evaluations, increases if i is added to Λ (work contribution) Then estimate the profit of each i as P(i) = ∆E(i) ∆W (i) and build the sparse grid using the set Λ of the M indices with the largest profit.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 28

slide-42
SLIDE 42

Stochastic multivariate polynomial approximation Optimized algorithms

Estimates for ∆W and ∆E

1

∆W (i): for nested points (e.g. Clenshaw Curtis, Gauss-Patterson) ∆W (i) = nb. new pts. in

N

  • n=1

∆m(in)

n

=

N

  • n=1

( m(in) − m(in − 1) )

2

∆E(i): we use the heuristic argument: use expansion on

  • rthnormal basis u =

p upψp

∆E(i) = ∆m(i)[u]V ⊗L2

ρ =

  • p

up∆m(i)ψpV ⊗L2

ρ

  • p≥m(i−1)

upV ∆m(i)ψpL2

ρ≈ um(i−1)V ∆m(i)ψm(i−1)L2 ρ

where um(i−1)V : estimated with a-priori / a-posteriori L(m(i)) := ∆m(i)ψm(i−1)L2

ρ: estimated numerically Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 29

slide-43
SLIDE 43

Stochastic multivariate polynomial approximation Optimized algorithms

Estimates for ∆W and ∆E

1

∆W (i): for nested points (e.g. Clenshaw Curtis, Gauss-Patterson) ∆W (i) = nb. new pts. in

N

  • n=1

∆m(in)

n

=

N

  • n=1

( m(in) − m(in − 1) )

2

∆E(i): we use the heuristic argument: use expansion on

  • rthnormal basis u =

p upψp

∆E(i) = ∆m(i)[u]V ⊗L2

ρ =

  • p

up∆m(i)ψpV ⊗L2

ρ

  • p≥m(i−1)

upV ∆m(i)ψpL2

ρ≈ um(i−1)V ∆m(i)ψm(i−1)L2 ρ

where um(i−1)V : estimated with a-priori / a-posteriori L(m(i)) := ∆m(i)ψm(i−1)L2

ρ: estimated numerically Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 29

slide-44
SLIDE 44

Stochastic multivariate polynomial approximation Optimized algorithms

Estimates for ∆W and ∆E

1

∆W (i): for nested points (e.g. Clenshaw Curtis, Gauss-Patterson) ∆W (i) = nb. new pts. in

N

  • n=1

∆m(in)

n

=

N

  • n=1

( m(in) − m(in − 1) )

2

∆E(i): we use the heuristic argument: use expansion on

  • rthnormal basis u =

p upψp

∆E(i) = ∆m(i)[u]V ⊗L2

ρ =

  • p

up∆m(i)ψpV ⊗L2

ρ

  • p≥m(i−1)

upV ∆m(i)ψpL2

ρ≈ um(i−1)V ∆m(i)ψm(i−1)L2 ρ

where um(i−1)V : estimated with a-priori / a-posteriori L(m(i)) := ∆m(i)ψm(i−1)L2

ρ: estimated numerically Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 29

slide-45
SLIDE 45

Stochastic multivariate polynomial approximation Optimized algorithms

All the pieces together

Optimal index set

Λ(w) =

  • i ∈ NN

+ : N

  • i=n

m(in − 1)gn − log |m(i − 1)|! m(i − 1)! −

N

  • n=1

log L(m(in)) m(in) − m(in − 1) ≤ w

  • (EW - Error Work grids)

where Legendre coeff + L(m(i)) = error estimate work estimate

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 30

slide-46
SLIDE 46

Stochastic multivariate polynomial approximation Optimized algorithms

A numerical check: a(x, y) = 1 + 0.1xy1 + 0.5x2y2

5 10 15 20 25 30 10

−15

10

−12

10

−9

10

−6

10

−3

10 ∆E(i) Legm(i−1)[ψ(u)] Legm(i−1)[ψ(u)] ⋅ Leb[m(i)]

Figure: Estimates of ∆E(i) for Θ(u)] in lexicographic order

20 40 60 80 100 120 140 10

−9

10

−7

10

−5

10

−3

10

−1

iso SM EW adaptive best M terms

Figure: Convergence plot for Θ(u) − Θ(uM)2

L2

ρ(Γ) w.r.t. |Λ|

Nested Clenshaw-Curtis knots The terms ∆E(i) have been computed with a sufficiently high level sparse grids. Comparison with dimension adaptive algorithm [Gerstner-Griebel ’03, Klimke, PhD ’06]

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 31

slide-47
SLIDE 47

Stochastic multivariate polynomial approximation Numerical results

Numerical test - 1D stationary lognormal field

L = 1, D = [0, L]2.      −∇ · a(y, x)∇u(y, x) = 0 u = 1 on x = 0, h = 0 on x = 1 no flux otherwise a(x, y) = eγ(x,y) µγ(x) = 0 Cov γ(x, x′) = σ2e−

|x1−x′ 1|2 LC2

We approximate γ as γ(y, x) ≈ µ(x) + σa0y0 + σ

K

  • k=1

ak

  • y2k−1 cos

π L kx1

  • + y2k sin

π L kx1

  • with yi ∼ N(0, 1), i.i.d.

Given the Fourier series σ2e− |z|2

LC2 = ∞

k=0 ck cos

π

L kz

  • , ak = √ck.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 32

slide-48
SLIDE 48

Stochastic multivariate polynomial approximation Numerical results

Numerical test - 1D stationary lognormal field

Quantity of interest: effective permeability E[Φ(u)], with Φ = L k(·, x)∂u(·, x) ∂x dx

  • Convergence: |E[Φ(uSG)] − E[Φ(u)]|

We compare Monte Carlo estimate with Knapsack grids based on Gauss-Hermite-Patterson points (nested Gauss-Hermite) Estimate of Hermite coefficients decay:

for the simpler problem ∇ · a(y)∇u(y, x) = f , a(y) = eb0+N

n=1 ynbn, we have uiV = C bin n

√in!.

Heuristic: use the same ansaz uiV ≈ C N

n=1 e−gnin √in! but

estimate the rates gn numerically.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 33

slide-49
SLIDE 49

Stochastic multivariate polynomial approximation Numerical results

Numerical test - 1D stationary lognormal field

Quantity of interest: effective permeability E[Φ(u)], with Φ = L k(·, x)∂u(·, x) ∂x dx

  • Convergence: |E[Φ(uSG)] − E[Φ(u)]|

We compare Monte Carlo estimate with Knapsack grids based on Gauss-Hermite-Patterson points (nested Gauss-Hermite) Estimate of Hermite coefficients decay:

for the simpler problem ∇ · a(y)∇u(y, x) = f , a(y) = eb0+N

n=1 ynbn, we have uiV = C bin n

√in!.

Heuristic: use the same ansaz uiV ≈ C N

n=1 e−gnin √in! but

estimate the rates gn numerically.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 33

slide-50
SLIDE 50

Stochastic multivariate polynomial approximation Numerical results

Numerical test - 1D stationary lognormal field

Quantity of interest: effective permeability E[Φ(u)], with Φ = L k(·, x)∂u(·, x) ∂x dx

  • Convergence: |E[Φ(uSG)] − E[Φ(u)]|

We compare Monte Carlo estimate with Knapsack grids based on Gauss-Hermite-Patterson points (nested Gauss-Hermite) Estimate of Hermite coefficients decay:

for the simpler problem ∇ · a(y)∇u(y, x) = f , a(y) = eb0+N

n=1 ynbn, we have uiV = C bin n

√in!.

Heuristic: use the same ansaz uiV ≈ C N

n=1 e−gnin √in! but

estimate the rates gn numerically.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 33

slide-51
SLIDE 51

Stochastic multivariate polynomial approximation Numerical results

Numerical test - 1D stationary lognormal field

Correlation length: LC = 0.2 Standard deviation: σ = 0.3 (c.o.v. ∼ 30%)

10 10

1

10

2

10

3

10

4

10

5

10

−10

10

−8

10

−6

10

−4

10

−2

10

4 11 14 17 19 21 24 29

sparse grid 1/N0.5 1/N1.7 MC run1 MC run2 MC run3 MC run4

The optimal set construction automatically adds new variables when needed. No need to truncate a-priori the random field

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 34

slide-52
SLIDE 52

Stochastic multivariate polynomial approximation Numerical results

How the sparse grid looks like

−5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 −5 5 η4 η3 η1 η2 η5 η6 η2 η3 η4 η5 η6 η7 Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 35

slide-53
SLIDE 53

Polynomial approximation by discrete projection on random points

Outline

1

Elliptic PDE with random coefficients

2

Stochastic multivariate polynomial approximation Galerkin projection Collocation on sparse grids Optimized algorithms Numerical results

3

Polynomial approximation by discrete projection on random points

4

Conclusions

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 38

slide-54
SLIDE 54

Polynomial approximation by discrete projection on random points

Discrete L2 projection using random evaluations (see poster 11 – G. Migliorati)

Another way to compute a polynomial approximation (besides Galerkin and sparse grid collocation) consists in doing a discrete least square approx. using random evaluations (see e.g. [Burkardt-Eldred

2009, Hosder-Walters et al. 2010, Eldred 2011, Blatman-Sudret 2008, ...])

1

Generate M random i.i.d. samples yk ∈ Γ, k = 1, . . . , M

2

Compute the corresponding solutions uk = u(y(k))

3

Find the discrete least square approximation ΠΛ,ω

M u ∈ V ⊗ PΛ(Γ)

ΠΛ,ω

M u = argmin v∈V ⊗PΛ(Γ)

1 M

M

  • k=1

uk − v(y(k))2

V

Two relevant questions For a given set Λ, how many samples should one use? What is the accuracy of the random discrete least square approx.?

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 39

slide-55
SLIDE 55

Polynomial approximation by discrete projection on random points

Discrete L2 projection using random evaluations (see poster 11 – G. Migliorati)

Another way to compute a polynomial approximation (besides Galerkin and sparse grid collocation) consists in doing a discrete least square approx. using random evaluations (see e.g. [Burkardt-Eldred

2009, Hosder-Walters et al. 2010, Eldred 2011, Blatman-Sudret 2008, ...])

1

Generate M random i.i.d. samples yk ∈ Γ, k = 1, . . . , M

2

Compute the corresponding solutions uk = u(y(k))

3

Find the discrete least square approximation ΠΛ,ω

M u ∈ V ⊗ PΛ(Γ)

ΠΛ,ω

M u = argmin v∈V ⊗PΛ(Γ)

1 M

M

  • k=1

uk − v(y(k))2

V

Two relevant questions For a given set Λ, how many samples should one use? What is the accuracy of the random discrete least square approx.?

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 39

slide-56
SLIDE 56

Polynomial approximation by discrete projection on random points

Discrete L2 projection using random evaluations (see poster 11 – G. Migliorati)

Another way to compute a polynomial approximation (besides Galerkin and sparse grid collocation) consists in doing a discrete least square approx. using random evaluations (see e.g. [Burkardt-Eldred

2009, Hosder-Walters et al. 2010, Eldred 2011, Blatman-Sudret 2008, ...])

1

Generate M random i.i.d. samples yk ∈ Γ, k = 1, . . . , M

2

Compute the corresponding solutions uk = u(y(k))

3

Find the discrete least square approximation ΠΛ,ω

M u ∈ V ⊗ PΛ(Γ)

ΠΛ,ω

M u = argmin v∈V ⊗PΛ(Γ)

1 M

M

  • k=1

uk − v(y(k))2

V

Two relevant questions For a given set Λ, how many samples should one use? What is the accuracy of the random discrete least square approx.?

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 39

slide-57
SLIDE 57

Polynomial approximation by discrete projection on random points

Some theoretical results [Migliorati-N.-von Schwerin-Tempone ’11]

For functions φ : Γ ⊂ RN → R, define continuous norm: φ2

L2

ρ =

  • Γ v 2(y)ρ(y)dy

discrete norm: φ2

M,ω = 1 M

M

i=1 φ(yi)2, with yi ∼ ρ(y)dy, i.i.d.

random discrete least square projection: ΠΛ,ω

M φ ∈ PΛ(Γ),

ΠΛ,ω

M φ = argmin v∈PΛ(Γ)

φ − v2

M,ω.

Theorem 1 Let C ω(M, Λ) := sup

v∈PΛ(Γ)

v2

L2

ρ

v2

M,ω

. Then

1

C ω(M, Λ) → 1 almost surely when M → ∞

2

φ − ΠΛ,ω

M φL2

ρ ≤ (1 +

  • C ω(M, Λ)) infv∈PΛ(Γ) φ − vL∞

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 40

slide-58
SLIDE 58

Polynomial approximation by discrete projection on random points

Some theoretical results [Migliorati-N.-von Schwerin-Tempone ’11]

So far we have only a result for 1D functions, Γ = [−1, 1], uniform distribution and approx. in the polynomial space Pw of degree w. Theorem 2 For any α ∈ (0, 1), let M be such that

M 3 log((M+1)/α) = 4

√ 3 w2 Then, it holds

P

  • φ − Πw,ω

M φL2

ρ ≤

  • 1 + 2
  • 3 log M + 1

α

  • inf

v∈Pw φ − vL∞

  • ≥ 1 − α.

Notice that log((M + 1)/α) ≈ log(Cw2/α) These results are confirmed numerically; the high dimentional case seems more foregiving w.r.t. the constraint M ∼ #Λ2 To achieve an approximation in PΛ, does this technique need less sampling points than the corresponding sparse grid ?

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 41

slide-59
SLIDE 59

Polynomial approximation by discrete projection on random points

Some theoretical results [Migliorati-N.-von Schwerin-Tempone ’11]

So far we have only a result for 1D functions, Γ = [−1, 1], uniform distribution and approx. in the polynomial space Pw of degree w. Theorem 2 For any α ∈ (0, 1), let M be such that

M 3 log((M+1)/α) = 4

√ 3 w2 Then, it holds

P

  • φ − Πw,ω

M φL2

ρ ≤

  • 1 + 2
  • 3 log M + 1

α

  • inf

v∈Pw φ − vL∞

  • ≥ 1 − α.

Notice that log((M + 1)/α) ≈ log(Cw2/α) These results are confirmed numerically; the high dimentional case seems more foregiving w.r.t. the constraint M ∼ #Λ2 To achieve an approximation in PΛ, does this technique need less sampling points than the corresponding sparse grid ?

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 41

slide-60
SLIDE 60

Polynomial approximation by discrete projection on random points

Some theoretical results [Migliorati-N.-von Schwerin-Tempone ’11]

So far we have only a result for 1D functions, Γ = [−1, 1], uniform distribution and approx. in the polynomial space Pw of degree w. Theorem 2 For any α ∈ (0, 1), let M be such that

M 3 log((M+1)/α) = 4

√ 3 w2 Then, it holds

P

  • φ − Πw,ω

M φL2

ρ ≤

  • 1 + 2
  • 3 log M + 1

α

  • inf

v∈Pw φ − vL∞

  • ≥ 1 − α.

Notice that log((M + 1)/α) ≈ log(Cw2/α) These results are confirmed numerically; the high dimentional case seems more foregiving w.r.t. the constraint M ∼ #Λ2 To achieve an approximation in PΛ, does this technique need less sampling points than the corresponding sparse grid ?

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 41

slide-61
SLIDE 61

Conclusions

Outline

1

Elliptic PDE with random coefficients

2

Stochastic multivariate polynomial approximation Galerkin projection Collocation on sparse grids Optimized algorithms Numerical results

3

Polynomial approximation by discrete projection on random points

4

Conclusions

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 48

slide-62
SLIDE 62

Conclusions

Conclusions

Solutions to elliptic equations with random coefficients typically feature analytic dependence on the parameters. Polynomial approximations are very effective. Sharp a-priori / a-posteriori analysis of the decay of the expansion of the solution in polynomial chaos allows to construct optimized polynomial spaces / sparse grids that provide effective approximations also in the infinite dimensional case. Discrete least square projection using random evaluations is a possible alternative to Galerkin or Collocation approaches. However, a better understanding is needed on the stability of the projection and the correct number of samples to use.

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 49

slide-63
SLIDE 63

Conclusions

References

  • G. Migliorati and F. Nobile and E. von Schwrin and R. Tempone

Analysis of the discrete L2 projection on polynomial spaces with random evaluations,

  • submitted. Available as MOX-Report.
  • J. B¨

ack and F. Nobile and L. Tamellini and R. Tempone On the optimal polynomial approximation of stochastic PDEs by Galerkin and Collocation methods, MOX Report 23/2011, to appear in M3AS.

  • I. Babuˇ

ska, F. Nobile and R. Tempone. A stochastic collocation method for elliptic PDEs with random input data, SIAM Review, 52(2):317–355, 2010

  • F. Nobile and R. Tempone

Analysis and implementation issues for the numerical approximation of parabolic equations with random coefficients, IJNME, 80:979–1006, 2009

  • F. Nobile, R. Tempone and C. Webster

An anisotropic sparse grid stochastic collocation method for PDEs with random input data, SIAM J. Numer. Anal., 46(5):2411–2442, 2008

Fabio Nobile (EPFL & PoliMi) Polynomial approximation of SPDEs RICAM, Linz, 12-16/12/2011 50