Application of quasi-Monte Carlo methods to PDEs with random coefficients Frances Kuo f.kuo@unsw.edu.au University of New South Wales, Sydney, Australia Built upon joint works with Ivan Graham, Rob Scheichl (Bath), Dirk Nuyens (KU Leuven), Christoph Schwab (ETH Zurich), Ian Sloan, James Nichols, Josef Dick, Thong Le Gia (UNSW). Frances Kuo
Motivating example Uncertainty in groundwater flow eg. risk analysis of radwaste disposal or CO 2 sequestration q + a � ∇ p = f Darcy’s law in D ⊂ R d , d = 1 , 2 , 3 ∇ · q = 0 mass conservation law together with boundary conditions [Cliffe, et. al. (2000)] x x x Uncertainty in a ( x x, ω ) leads to uncertainty in q ( x x, ω ) and p ( x x, ω ) Frances Kuo
Expected values of quantities of interest To compute the expected value of some quantity of interest: 1. Generate a number of realizations of the random field (Some approximation may be required) 2. For each realization, solve the PDE using e.g. FEM / FVM / mFEM 3. Take the average of all solutions from different realizations This describes Monte Carlo simulation. Example : particle dispersion ∂p n = 0 ∂� 1 0.8 p = 1 p = 0 0.6 release point → ◦ 0.4 • 0.2 0 0 0.2 0.4 0.6 0.8 1 ∂p n = 0 ∂� Frances Kuo
Expected values of quantities of interest To compute the expected value of some quantity of interest: 1. Generate a number of realizations of the random field (Some approximation may be required) 2. For each realization, solve the PDE using e.g. FEM / FVM / mFEM 3. Take the average of all solutions from different realizations This describes Monte Carlo simulation. NOTE : expected value = (high dimensional) integral → use quasi-Monte Carlo methods −1 10 MC n − 1 / 2 −2 10 −3 10 error n − 1 or better −4 10 QMC s = stochastic dimension −5 10 3 4 5 6 10 10 10 10 n Frances Kuo
MC v.s. QMC in the unit cube n � 1 � [0 , 1] s F ( y y y ) d y y y ≈ F ( t t t i ) n i =1 Monte Carlo method Quasi-Monte Carlo methods t t t t i deterministic t t i random uniform close to n − 1 convergence or better n − 1 / 2 convergence more effective for earlier variables and lower -order projections order of variables irrelevant order of variables very important 1 1 1 • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • 0 1 0 1 0 1 First 64 points of a 64 random points A lattice rule with 64 points 2D Sobol ′ sequence use randomized QMC methods for error estimation Frances Kuo
QMC Two main families of QMC methods: (t,m,s)-nets and (t,s)-sequences Niederreiter book (1992) Sloan and Joe book (1994) lattice rules Important developments: component-by-component ( CBC ) A group under addition modulo Z construction and includes the integer points higher order digital nets • • • • • • • • • • • • • • • 1 1 • • • • • • • Dick and Pillichshammer book (2010) • • • • • • • • • • • • • • • • • • • • • Dick, Kuo, Sloan Acta Numerica (2013) • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • Having the right number of • 0 1 0 1 • • • First 64 points of a • • A lattice rule with 64 points points in various sub-cubes 2D Sobol ′ sequence • (0,6,2)-net Frances Kuo
Application of QMC theory Worst case error bound � � � � n n � � � � y − 1 y − 1 � � � � � ≤ e wor � ≤ e wor ( t � � � � y y y y t t t t t t [0 , 1] s F ( y [0 , 1] s F ( y y ) d y y ) d y F ( t F ( t t i ) t i ) ( t t 1 , . . . ,t t 1 , . . . ,t t n ) � F � γ t n ) � F � γ γ � � � � γ γ γ n n � � i =1 i =1 Weighted Sobolev space [Sloan & Wo´ zniakowski, 1998] � � 2 � ∂ | u | F 1 � � � � F � 2 � � y y γ = ( y y u ; 0) d y y u � � γ γ y γ u ∂y y u [0 , 1] | u | � � u ⊆{ 1: s } “anchor” at 0 (also “unanchored”) 2 s subsets “weights” Mixed first derivatives are square integrable Small weight γ u means that F depends weakly on the variables x x x u Choose weights to minimize the error bound � b u � � 1 / (2 λ ) � � 1 / 2 � 1 / (1+ λ ) 2 b u � � γ λ u a u ⇒ γ u = n γ u a u u ⊆{ 1: s } u ⊆{ 1: s } � �� � � �� � bound on worst case error (CBC) bound on norm Construct points (CBC) to minimize the worst case error Frances Kuo
The present state of play ... in the application of QMC to PDEs with random coefficients: 1. Graham, K., Nuyens, Scheichl, Sloan (J. Comput. Physics, 2011) Application of QMC to the lognormal case, no error analysis Use circulant embedding to avoid truncation of KL expansion 2. K., Schwab, Sloan (SIAM J. Numer. Anal., 2012) Application of QMC to the uniform case [cf. Cohen, De Vore, Schwab, 2010] First complete error analysis: how to choose “weights” 3. K., Schwab, Sloan (submitted) [cf. Heinrich 1998; Giles 2008] A multi-level version of the analysis for the uniform case 4. Graham, K., Nichols, Scheichl, Schwab, Sloan (submitted) Application of QMC to the lognormal case 5. Schwab (Proceedings of MCQMC 2012) Application of QMC to the uniform case for affine operator equations 6. Le Gia (Proceedings of MCQMC 2012) Application of QMC to the uniform case for PDE over the sphere Frances Kuo
The present state of play (continued) 7. Dick, K., Le Gia, Nuyens, Schwab (SIAM J. Numer. Anal., to appear) Application of higher order QMC to the uniform case for affine operator equations New Banach (non-Hilbert) space setting New SPOD weights(smoothness-driven product and order dependent) New fast CBC construction of interlaced polynomial lattice rules 8. Dick, K., Le Gia, Schwab (submitted) A multi-level version of the higher order analysis for the uniform case . . . Lots of related works: http://people.maths.ox.ac.uk/~gilesm/mlmc_community.html Also QMC works by Harbrecht, Peters, Siebenmorgen Frances Kuo
The uniform case [K., Schwab, Sloan (2012)] x y x y x x y −∇ · ( a ( x x,y y ) ∇ u ( x x,y y )) = f ( x x ) in D , u ( x x,y y ) = 0 on ∂D ∞ � y ∈ ( − 1 2 , 1 2 ) ∞ x y x x y a ( x x,y y ) = ¯ a ( x x ) + y j ψ j ( x x ) , y j =1 Assumptions: x y There exists a min and a max such that 0 < a min ≤ a ( x x,y y ) ≤ a max There exists p 0 ∈ (0 , 1] such that � ∞ j =1 � ψ j � p 0 L ∞ ( D ) < ∞ · · · y y Goal: for F ( y y ) = G ( u ( · ,y y )) a linear functional of u , we want to approximate � � y y 2 ) ∞ F ( y y ) d y y = lim 2 ) s F ( y 1 , . . . , y s , 0 , 0 , . . . ) d y 1 · · · d y s ( − 1 2 , 1 ( − 1 2 , 1 s →∞ Our strategy: Truncate the infinite sum to s terms (truncation error) Solve the PDE using FEM (discretization error) Approximate the integral using QMC (quadrature error) Frances Kuo
Error analysis – a sum of three terms � ∞ j =1 � ψ j � p 0 I ( G ( u )) ≈ Q s,n ( G ( u s h )) L ∞ ( D ) < ∞ n f ∈ H − 1+ t ( D ) = 1 � G ( u s t h ( · ,t t i )) G ∈ H − 1+ t ′ ( D ) n i =1 = O ( s − 2(1 /p 0 − 1) ) Truncation error use “integration” Discretization error = O ( h t + t ′ ) Aubin-Nitsche duality trick = O ( n − min(1 /p 0 − 1 / 2 , 1 − δ ) ) Quadrature error constant independent of s y Differentiate the PDE w.r.t. y y [Cohen, De Vore, Schwab, 2010] ν | ! � j ≥ 1 � ψ j � ν j ν � ∂ ν ν y ν y u ( · ,y y ) � H 1 0 ( D ) ≤ C | ν y L ∞ ( D ) y Weighted Sobolev space with square-integrable mixed first derivatives Fast CBC construction of randomly shifted lattice rules Minimize error bound: get POD weights (product and order dependent) � � b c | u | | u | ! � γ u = j ∈ u � ψ j � L ∞ ( D ) Frances Kuo
Recommend
More recommend
