Richardson extrapolation and higher order QMC Takashi Goda School - - PowerPoint PPT Presentation

▶
richardson extrapolation and higher order qmc
SMART_READER_LITE
LIVE PREVIEW

Richardson extrapolation and higher order QMC Takashi Goda School - - PowerPoint PPT Presentation

Jun 22, 2023 419 likes •1.83k views

Richardson extrapolation and higher order QMC Takashi Goda School of Engineering, University of Tokyo MCM 2019 Partly joint with Josef Dick and Takehito Yoshiki Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 1 / 61 Outline QMC

slide-1
SLIDE 1

Richardson extrapolation and higher order QMC

Takashi Goda

School of Engineering, University of Tokyo

MCM 2019 Partly joint with Josef Dick and Takehito Yoshiki

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 1 / 61

slide-2
SLIDE 2

Outline

QMC and digital nets/sequences Classical QMC Higher order QMC Richardson extrapolation and QMC Application 1: Truncation of higher order nets and sequences Application 2: Extrapolated polynomial lattice rules

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 2 / 61

slide-3
SLIDE 3

Outline

QMC and digital nets/sequences Classical QMC Higher order QMC Richardson extrapolation and QMC Application 1: Truncation of higher order nets and sequences Application 2: Extrapolated polynomial lattice rules

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 3 / 61

slide-4
SLIDE 4

Quasi-Monte Carlo

Approximate I(f ) = ∫

[0,1)s f (x) dx,

where f : [0, 1)s → R is integrable. Choose an N-element point set P ⊂ [0, 1)s and approximate I(f ) by QP(f ) = 1 N ∑

x∈P

f (x). This is called a quasi-Monte Carlo (QMC) integration over P.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 4 / 61

slide-5
SLIDE 5

Digital nets in base 2

Let m, n ∈ N := {1, 2, . . .}. Let C1, . . . , Cs ∈ Fn×m

2

(F2 := {0, 1}, the two-element field).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 5 / 61

slide-6
SLIDE 6

Digital nets in base 2

Let m, n ∈ N := {1, 2, . . .}. Let C1, . . . , Cs ∈ Fn×m

2

(F2 := {0, 1}, the two-element field). Denote the dyadic expansion of 0 ≤ h < 2m by h = (ηm−1 . . . η1η0)2.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 5 / 61

slide-7
SLIDE 7

Digital nets in base 2

Let m, n ∈ N := {1, 2, . . .}. Let C1, . . . , Cs ∈ Fn×m

2

(F2 := {0, 1}, the two-element field). Denote the dyadic expansion of 0 ≤ h < 2m by h = (ηm−1 . . . η1η0)2. For 1 ≤ j ≤ s, let xj = (0.ξ(j)

1 . . . ξ(j) n 00 . . .)2,

where    ξ(j)

1

. . . ξ(j)

n

   = Cj    η0 . . . ηm−1    ∈ Fn

2.

This gives a point x = (x1, . . . , xs) ∈ [0, 1)s.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 5 / 61

slide-8
SLIDE 8

Digital nets in base 2

Let m, n ∈ N := {1, 2, . . .}. Let C1, . . . , Cs ∈ Fn×m

2

(F2 := {0, 1}, the two-element field). Denote the dyadic expansion of 0 ≤ h < 2m by h = (ηm−1 . . . η1η0)2. For 1 ≤ j ≤ s, let xj = (0.ξ(j)

1 . . . ξ(j) n 00 . . .)2,

where    ξ(j)

1

. . . ξ(j)

n

   = Cj    η0 . . . ηm−1    ∈ Fn

2.

This gives a point x = (x1, . . . , xs) ∈ [0, 1)s. Running through all 0 ≤ h < 2m, we get a 2m-element point set P called a digital net over F2.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 5 / 61

slide-9
SLIDE 9

Digital sequences in base 2

Let C1, . . . , Cs ∈ FN×N

2

, where every column is assumed to have only finite non-zero elements. Denote the dyadic expansion of h ∈ N0 := N ∪ {0} by h = (. . . η1η0)2. where all but finite number of ηi are 0. For 1 ≤ j ≤ s, let xj = (0.ξ(j)

1 ξ(j) 2 . . .)2,

where    ξ(j)

1

ξ(j)

2

. . .    = Cj    η0 η1 . . .    ∈ FN

2 .

This gives a point x = (x1, . . . , xs) ∈ [0, 1)s. This way we get an infinite sequence of points S called a digital sequence over F2.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 6 / 61

slide-10
SLIDE 10

Some remarks

“Goodness” (undefined) of digital nets/sequences is determined by generating matrices C1, . . . , Cs. For nets, m determines total number of points, while n determines precision of each point. For sequences, the assumption that every column has only finite non-zero elements is to ensure that each point has finite precision.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 7 / 61

slide-11
SLIDE 11

Precision?

For classical QMC where “goodness” is measured by discrepancy, it is usual to set n = m (square matrices) for nets, or consider upper-triangular matrices for sequences. For higher order QMC where “goodness” is measured by the worst-case error for a class of smooth functions, we need higher

  • precision. Typically, n = αm for nets (α: dominating mixed

smoothness of functions).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 8 / 61

slide-12
SLIDE 12

Precision?

For classical QMC where “goodness” is measured by discrepancy, it is usual to set n = m (square matrices) for nets, or consider upper-triangular matrices for sequences. For higher order QMC where “goodness” is measured by the worst-case error for a class of smooth functions, we need higher

  • precision. Typically, n = αm for nets (α: dominating mixed

smoothness of functions).

Aim of this talk

is to show that Richardson extrapolation can be used as a technique to reduce necessary precision of higher order QMC from αm to m.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 8 / 61

slide-13
SLIDE 13

Precision?

For classical QMC where “goodness” is measured by discrepancy, it is usual to set n = m (square matrices) for nets, or consider upper-triangular matrices for sequences. For higher order QMC where “goodness” is measured by the worst-case error for a class of smooth functions, we need higher

  • precision. Typically, n = αm for nets (α: dominating mixed

smoothness of functions).

Aim of this talk

is to show that Richardson extrapolation can be used as a technique to reduce necessary precision of higher order QMC from αm to m.

Implication

Possibility to try large values of α w/o suffering from round-off error. Search space for C1, . . . , Cs is significantly reduced.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 8 / 61

slide-14
SLIDE 14

Outline

QMC and digital nets/sequences Classical QMC Higher order QMC Richardson extrapolation and QMC Application 1: Truncation of higher order nets and sequences Application 2: Extrapolated polynomial lattice rules

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 9 / 61

slide-15
SLIDE 15

Local discrepancy function

For an N-element point set P and y ∈ [0, 1)s, let ∆P(y) := 1 N ∑

x∈P

1x∈[0,y) − λ([0, y)), where λ denotes the Lebesgue measure.

y

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 10 / 61

slide-16
SLIDE 16

Discrepancy and Koksma-Hlawka inequality

For 1 ≤ p ≤ ∞, the Lp-discrepancy of P is defined by Lp(P) := (∫

[0,1)s |∆P(y)|p dy

)1/p . When p = ∞, the obvious modification is needed.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 11 / 61

slide-17
SLIDE 17

Discrepancy and Koksma-Hlawka inequality

For 1 ≤ p ≤ ∞, the Lp-discrepancy of P is defined by Lp(P) := (∫

[0,1)s |∆P(y)|p dy

)1/p . When p = ∞, the obvious modification is needed. The integration error is bounded by |QP(f ) − I(f )| ≤ VHK(f )L∞(P) where VHK(f ) is the Hardy-Krause variation of f .

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 11 / 61

slide-18
SLIDE 18

Discrepancy and Koksma-Hlawka inequality

For 1 ≤ p ≤ ∞, the Lp-discrepancy of P is defined by Lp(P) := (∫

[0,1)s |∆P(y)|p dy

)1/p . When p = ∞, the obvious modification is needed. The integration error is bounded by |QP(f ) − I(f )| ≤ VHK(f )L∞(P) where VHK(f ) is the Hardy-Krause variation of f .

Main focus on classical QMC

How to construct digital nets/sequences with small L∞(P).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 11 / 61

slide-19
SLIDE 19

Quality measure: t-value

Let t be an integer such that, for any choice d1, . . . , ds ≥ 0 with d1 + · · · + ds = m − t, the first d1 rows of C1 . . . the first ds rows of Cs are linearly independent over F2. P is called a digital (t, m, s)-net.

... C1 C2 Cs

d1 d2 ds

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 12 / 61

slide-20
SLIDE 20

Quality measure: t-value

An equi-distribution property of digital (t, m, s)-nets: every dyadic elementary box of the form

s

∏

j=1

[ aj 2cj , aj + 1 2cj ) with      c1, . . . , cs ≥ 0, c1 + · · · + cs = m − t, 0 ≤ aj < 2cj contains exactly 2t points.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 13 / 61

slide-21
SLIDE 21

Quality measure: t-value

An equi-distribution property of digital (t, m, s)-nets: every dyadic elementary box of the form

s

∏

j=1

[ aj 2cj , aj + 1 2cj ) with      c1, . . . , cs ≥ 0, c1 + · · · + cs = m − t, 0 ≤ aj < 2cj contains exactly 2t points. Let t be an integer such that, for any m ≥ t, the first 2m points of S are a digital (t, m, s)-net. S is called a digital (t, s)-sequence.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 13 / 61

slide-22
SLIDE 22

t-value and star-discrepancy

If P is a digital (t, m, s)-net over F2, L∞(P) ≤ Cs,t (log N)s−1 N with N = 2m. If P is the first N points of a digital (t, s)-sequence over F2, L∞(P) ≤ Ds,t (log N)s N , for all N ≥ 2.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 14 / 61

slide-23
SLIDE 23

Explicit construction

For s = 1, C1 can be the identity matrix, which generates the famous digital (0, 1)-sequence called van der Corput sequence.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 15 / 61

slide-24
SLIDE 24

Explicit construction

For s = 1, C1 can be the identity matrix, which generates the famous digital (0, 1)-sequence called van der Corput sequence. For s = 2, C1 and C2 can be C1 =      1 · · · 1 · · · . . . . . . ... . . . · · · 1      , C2 =      · · · 1 · · · 1 . . . ... . . . . . . 1 · · ·      , which generates a digital (0, m, 2)-net, known as the Hammersley point set. Not extensible in m.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 15 / 61

slide-25
SLIDE 25

Explicit construction

Figure: Hammersley point set for m = 6

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 16 / 61

slide-26
SLIDE 26

Explicit construction

Figure: Hammersley point set for m = 6

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 16 / 61

slide-27
SLIDE 27

Explicit construction

Figure: Hammersley point set for m = 6

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 16 / 61

slide-28
SLIDE 28

Explicit construction

Figure: Hammersley point set for m = 6

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 16 / 61

slide-29
SLIDE 29

Explicit construction

Figure: Hammersley point set for m = 6

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 16 / 61

slide-30
SLIDE 30

Explicit construction

Figure: Hammersley point set for m = 6

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 16 / 61

slide-31
SLIDE 31

Explicit construction

Figure: Hammersley point set for m = 6

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 16 / 61

slide-32
SLIDE 32

Explicit construction

Figure: Hammersley point set for m = 6

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 16 / 61

slide-33
SLIDE 33

Explicit construction

Let p1, p2, . . . ∈ F2[x] be a sequence of distinct primitive/irreducible polynomials over F2 with e1 ≤ e2 ≤ · · · where ej = deg(pj).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 17 / 61

slide-34
SLIDE 34

Explicit construction

Let p1, p2, . . . ∈ F2[x] be a sequence of distinct primitive/irreducible polynomials over F2 with e1 ≤ e2 ≤ · · · where ej = deg(pj). For each j, Cj = (c(j)

k,l ) is given by the coefficients of the following

Laurent series: xej−1 pj(x) = c(j)

1,1

x + c(j)

1,2

x2 + · · · . . . 1 pj(x) = c(j)

ej,1

x + c(j)

ej,2

x2 + · · · xej−1 (pj(x))2 = c(j)

ej+1,1

x + c(j)

ej+1,2

x2 + · · · . . . (Sobol’ 1967; Niederreiter 1988; Tezuka 1993)

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 17 / 61

slide-35
SLIDE 35

Explicit construction

All matrices Cj ∈ FN×N

2

are upper triangular and generate a digital (t, s)-sequence with t = (e1 − 1) + · · · + (es − 1). (I used a C implementation for this sequence in at least more than 58636 dimensions due to Tomohiko Hironaka.) Figure: 2D projections of the first 26 points of the Niederreiter sequence.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 18 / 61

slide-36
SLIDE 36

Polynomial lattice point sets (Niederreiter, 1992)

Let p ∈ F2[x] be irreducible with deg(p) = m Let q = (q1, . . . , qs) ∈ (F2[x])s with deg(qj) < m.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 19 / 61

slide-37
SLIDE 37

Polynomial lattice point sets (Niederreiter, 1992)

Let p ∈ F2[x] be irreducible with deg(p) = m Let q = (q1, . . . , qs) ∈ (F2[x])s with deg(qj) < m. For each j, Cj is given by the square Hankel matrix Cj =       a(j)

1

a(j)

2

· · · a(j)

m

a(j)

2

a(j)

3

. . . a(j)

m+1

. . . . . . ... . . . a(j)

m

a(j)

m+1

. . . a(j)

2m−1

      ∈ Fm×m

2

where a(j)

1 , a(j) 2 , . . . are the coefficients of the Laurent series

qj(x) p(x) = a(j)

1

x + a(j)

2

x2 + · · · .

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 19 / 61

slide-38
SLIDE 38

Polynomial lattice point sets (Niederreiter, 1992)

Let p ∈ F2[x] be irreducible with deg(p) = m Let q = (q1, . . . , qs) ∈ (F2[x])s with deg(qj) < m. For each j, Cj is given by the square Hankel matrix Cj =       a(j)

1

a(j)

2

· · · a(j)

m

a(j)

2

a(j)

3

. . . a(j)

m+1

. . . . . . ... . . . a(j)

m

a(j)

m+1

. . . a(j)

2m−1

      ∈ Fm×m

2

where a(j)

1 , a(j) 2 , . . . are the coefficients of the Laurent series

qj(x) p(x) = a(j)

1

x + a(j)

2

x2 + · · · . The resulting digital net is called a polynomial lattice point set P(p, q). Usually the vector q is constructed by a (fast) computer search algorithm (Nuyens & Cools, 2006;

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 19 / 61

slide-39
SLIDE 39

Polynomial lattice point sets (Niederreiter, 1992)

Let p ∈ F2[x] be irreducible with deg(p) = m Let q = (q1, . . . , qs) ∈ (F2[x])s with deg(qj) < m. For each j, Cj is given by the square Hankel matrix Cj =       a(j)

1

a(j)

2

· · · a(j)

m

a(j)

2

a(j)

3

. . . a(j)

m+1

. . . . . . ... . . . a(j)

m

a(j)

m+1

. . . a(j)

2m−1

      ∈ Fm×m

2

where a(j)

1 , a(j) 2 , . . . are the coefficients of the Laurent series

qj(x) p(x) = a(j)

1

x + a(j)

2

x2 + · · · . The resulting digital net is called a polynomial lattice point set P(p, q). Usually the vector q is constructed by a (fast) computer search algorithm (Nuyens & Cools, 2006; P. Kritzer, this afternoon!).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 19 / 61

slide-40
SLIDE 40

Outline

QMC and digital nets/sequences Classical QMC Higher order QMC Richardson extrapolation and QMC Application 1: Truncation of higher order nets and sequences Application 2: Extrapolated polynomial lattice rules

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 20 / 61

slide-41
SLIDE 41

What if f is smooth

In some applications, such as PDEs with random coefficients, f can be smooth. A proper design of QMC point sets enables higher order convergence of the integration error than O(1/N) as expected from the KH inequality. So far, QMC point sets achieving higher order convergence for non-periodic smooth functions are

1 Higher order digital nets/sequences (Dick, 2008; ...), and 2 Tent-transformed lattice point sets (Hickernell, 2002; ...). Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 21 / 61

slide-42
SLIDE 42

What if f is smooth

In some applications, such as PDEs with random coefficients, f can be smooth. A proper design of QMC point sets enables higher order convergence of the integration error than O(1/N) as expected from the KH inequality. So far, QMC point sets achieving higher order convergence for non-periodic smooth functions are

1 Higher order digital nets/sequences (Dick, 2008; ...), and 2 Tent-transformed lattice point sets (Hickernell, 2002; ...).

In this talk, I will focus on higher order digital nets/sequences. The contents of this section are mostly developed by Dick (2008). Please refer to a recent review by G. & Suzuki (arXiv:1903.12353).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 21 / 61

slide-43
SLIDE 43

Quality measure: high order t-value

Let α ∈ N. Let t be an integer such that, for any choice 1 ≤ dj,vj < · · · < dj,1 ≤ αm, 0 ≤ vj ≤ αm, 1 ≤ j ≤ s, with

s

∑

j=1 min(vj,α)

∑

i=1

dj,i = αm − t, the d1,v1, . . . , d1,1-th rows of C1 . . . the ds,vs, . . . , ds,1-th rows of Cs are linearly independent over F2. P is called an order α digital (t, m, s)-net.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 22 / 61

slide-44
SLIDE 44

Quality measure: high order t-value

Order α digital (t, m, s)-nets hold equi-distribution properties: union

  • f dyadic elementary boxes contains the fair number of points (shown

later visually).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 23 / 61

slide-45
SLIDE 45

Quality measure: high order t-value

Order α digital (t, m, s)-nets hold equi-distribution properties: union

  • f dyadic elementary boxes contains the fair number of points (shown

later visually). Let t be an integer such that, for any αm ≥ t, the first 2m points of S are an order α digital (t, m, s)-net. S is called an order α digital (t, s)-sequence.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 23 / 61

slide-46
SLIDE 46

Explicit construction

Define Dα : [0, 1)α → [0, 1) by            x1 = (0.ξ(1)

1 ξ(1) 2 ξ(1) 3

. . .)2 x2 = (0.ξ(2)

1 ξ(2) 2 ξ(2) 3

. . .)2 . . . xα = (0.ξ(α)

1

ξ(α)

2

ξ(α)

3

. . .)2 → (0. ξ(1)

1 ξ(2) 1

. . . ξ(α)

1

  • α

ξ(1)

2 ξ(2) 2

. . . ξ(α)

2

  • α

. . .)2.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 24 / 61

slide-47
SLIDE 47

Explicit construction

Define Dα : [0, 1)α → [0, 1) by            x1 = (0.ξ(1)

1 ξ(1) 2 ξ(1) 3

. . .)2 x2 = (0.ξ(2)

1 ξ(2) 2 ξ(2) 3

. . .)2 . . . xα = (0.ξ(α)

1

ξ(α)

2

ξ(α)

3

. . .)2 → (0. ξ(1)

1 ξ(2) 1

. . . ξ(α)

1

  • α

ξ(1)

2 ξ(2) 2

. . . ξ(α)

2

  • α

. . .)2. For x ∈ [0, 1)αs, let Dα(x) = (Dα(x1, . . . , xα), Dα(xα+1, . . . , x2α), . . . , Dα(xα(s−1)+1, . . . , xαs)) ∈ [0, 1)s.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 24 / 61

slide-48
SLIDE 48

Explicit construction

Define Dα : [0, 1)α → [0, 1) by            x1 = (0.ξ(1)

1 ξ(1) 2 ξ(1) 3

. . .)2 x2 = (0.ξ(2)

1 ξ(2) 2 ξ(2) 3

. . .)2 . . . xα = (0.ξ(α)

1

ξ(α)

2

ξ(α)

3

. . .)2 → (0. ξ(1)

1 ξ(2) 1

. . . ξ(α)

1

  • α

ξ(1)

2 ξ(2) 2

. . . ξ(α)

2

  • α

. . .)2. For x ∈ [0, 1)αs, let Dα(x) = (Dα(x1, . . . , xα), Dα(xα+1, . . . , x2α), . . . , Dα(xα(s−1)+1, . . . , xαs)) ∈ [0, 1)s. For a digital (t, m, αs)-net P, we write Dα(P) = {Dα(x) | x ∈ P}.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 24 / 61

slide-49
SLIDE 49

Another look at construction

Let P be a digital (t, m, αs)-net with C1, . . . , Cαs ∈ Fm×m

2

. We write

C1 =    c(1)

1

. . . c(1)

m

   , C2 =    c(2)

1

. . . c(2)

m

   , . . . , Cα =    c(α)

1

. . . c(α)

m

   , . . . .

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 25 / 61

slide-50
SLIDE 50

Another look at construction

Let P be a digital (t, m, αs)-net with C1, . . . , Cαs ∈ Fm×m

2

. We write

C1 =    c(1)

1

. . . c(1)

m

   , C2 =    c(2)

1

. . . c(2)

m

   , . . . , Cα =    c(α)

1

. . . c(α)

m

   , . . . .

Dα(P) is a digital net with D1, . . . , Ds ∈ Fαm×m

2

where

D1 =                    c(1)

1

c(2)

1

. . . c(α)

1

. . . c(1)

m

c(2)

m

. . . c(α)

m

                   , D2 =                    c(α+1)

1

c(α+2)

1 .

. . c(2α)

1.

. . c(α+1)

m

c(α+2)

m

. . . c(2α)

m

                   , . . . .

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 25 / 61

slide-51
SLIDE 51

High order t-value

For a digital (t, m, αs)-net P, Dα(P) is an order α digital (t′, m, s)-net with t′ ≤ α min { m, t + ⌊s(α − 1) 2 ⌋} . For a digital (t, αs)-sequences S, Dα(S) is an order α digital (t′, s)-sequences with t′ ≤ αt + sα(α − 1) 2 .

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 26 / 61

slide-52
SLIDE 52

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-53
SLIDE 53

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-54
SLIDE 54

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-55
SLIDE 55

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-56
SLIDE 56

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-57
SLIDE 57

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-58
SLIDE 58

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-59
SLIDE 59

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-60
SLIDE 60

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-61
SLIDE 61

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-62
SLIDE 62

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-63
SLIDE 63

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-64
SLIDE 64

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-65
SLIDE 65

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-66
SLIDE 66

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-67
SLIDE 67

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-68
SLIDE 68

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-69
SLIDE 69

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-70
SLIDE 70

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-71
SLIDE 71

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-72
SLIDE 72

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-73
SLIDE 73

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-74
SLIDE 74

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-75
SLIDE 75

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-76
SLIDE 76

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-77
SLIDE 77

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-78
SLIDE 78

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-79
SLIDE 79

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-80
SLIDE 80

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-81
SLIDE 81

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-82
SLIDE 82

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-83
SLIDE 83

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-84
SLIDE 84

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-85
SLIDE 85

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-86
SLIDE 86

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-87
SLIDE 87

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-88
SLIDE 88

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-89
SLIDE 89

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-90
SLIDE 90

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-91
SLIDE 91

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-92
SLIDE 92

Interlaced Sobol’ point sets with α = 2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 27 / 61

slide-93
SLIDE 93

Interlaced polynomial lattice point sets

Let p ∈ F2[x] be irreducible with deg(p) = m, and let q = (q1, . . . , qαs) ∈ (F2[x])αs with deg(qj) < m. An interlaced polynomial lattice point set is just Dα(P(p, q)).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 28 / 61

slide-94
SLIDE 94

Interlaced polynomial lattice point sets

Let p ∈ F2[x] be irreducible with deg(p) = m, and let q = (q1, . . . , qαs) ∈ (F2[x])αs with deg(qj) < m. An interlaced polynomial lattice point set is just Dα(P(p, q)). The vector q can be constructed component by component. In the simplest case, the necessary construction cost is of O(αsN log N) with O(N) memory (G. & Dick, 2015; G., 2015). In applications to PDEs with random coefficients, the criterion sometimes becomes a bit complicated, requiring O(αsN log N + α2s2N) construction cost with O(αsN) memory (Dick, Kuo, Le Gia, Nuyens & Schwab, 2014).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 28 / 61

slide-95
SLIDE 95

Outline

QMC and digital nets/sequences Classical QMC Higher order QMC Richardson extrapolation and QMC Application 1: Truncation of higher order nets and sequences Application 2: Extrapolated polynomial lattice rules

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 29 / 61

slide-96
SLIDE 96

Walsh functions

For k = (. . . κ1κ0)2 ∈ N0, the k-th Walsh function is defined by walk(x) = (−1)κ0ξ1+κ1ξ2+···, where x = (0.ξ1ξ2 . . .)2 ∈ [0, 1), unique in the sense that infinitely many ξi are equal to 0.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 30 / 61

slide-97
SLIDE 97

Walsh functions

For k = (. . . κ1κ0)2 ∈ N0, the k-th Walsh function is defined by walk(x) = (−1)κ0ξ1+κ1ξ2+···, where x = (0.ξ1ξ2 . . .)2 ∈ [0, 1), unique in the sense that infinitely many ξi are equal to 0. For s ≥ 1 and k = (k1, . . . , ks) ∈ Ns

0, the k-th Walsh function is

defined by walk(x) =

s

∏

j=1

walkj(xj).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 30 / 61

slide-98
SLIDE 98

Walsh functions

Figure: The k-th Walsh functions for k = 0, 1, 2, 3

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 31 / 61

slide-99
SLIDE 99

Walsh functions

Figure: The k-th Walsh functions for k = 4, 5, 6, 7

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 32 / 61

slide-100
SLIDE 100

Walsh functions

Every Walsh function is a piecewise constant function.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 33 / 61

slide-101
SLIDE 101

Walsh functions

Every Walsh function is a piecewise constant function. Importantly the system of Walsh functions {walk | k ∈ Ns

0} is a

complete orthonormal basis in L2([0, 1)s). Thus we have a Walsh series of f : ∑

k∈Ns

ˆ f (k)walk(x) where ˆ f (k) = ∫

[0,1)s f (x)walk(x) dx.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 33 / 61

slide-102
SLIDE 102

Dual nets

For a digital net P with C1, . . . , Cs ∈ Fn×m

2

, the dual net is defined by P⊥ =      k ∈ Ns

  • C ⊤

1

   κ1,0 . . . κ1,n−1    ⊕ · · · ⊕ C ⊤

s

   κs,0 . . . κs,n−1    = 0 ∈ Fm

2

     , where we write kj = (. . . κj,1κj,0)2.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 34 / 61

slide-103
SLIDE 103

Dual nets

For a digital net P with C1, . . . , Cs ∈ Fn×m

2

, the dual net is defined by P⊥ =      k ∈ Ns

  • C ⊤

1

   κ1,0 . . . κ1,n−1    ⊕ · · · ⊕ C ⊤

s

   κs,0 . . . κs,n−1    = 0 ∈ Fm

2

     , where we write kj = (. . . κj,1κj,0)2.

Trivial fact

k ∈ P⊥ if 2n | k, i.e., 2n | kj for all j = 1, . . . , s.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 34 / 61

slide-104
SLIDE 104

Character property

Most of the Walsh functions can be integrated exactly: 1 2m ∑

x∈P

walk(x) = { 1 if k ∈ P⊥,

  • therwise.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 35 / 61

slide-105
SLIDE 105

Character property

Most of the Walsh functions can be integrated exactly: 1 2m ∑

x∈P

walk(x) = { 1 if k ∈ P⊥,

  • therwise.

I will also use the following fact: 1 2n

2n−1

∑

i=0

walk ( i 2n ) = 1 2n ∑

ξ0,...,ξn−1∈{0,1} n

∏

i=1

(−1)κi−1ξi = 1 2n

n

∏

i=1

∑

ξi∈{0,1}

(−1)κi−1ξi = { 1 if κ0 = . . . = κn−1 = 0, i.e., if 2n | k,

  • therwise.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 35 / 61

slide-106
SLIDE 106

Decomposition of integration error

For a digital net P with C1, . . . , Cs ∈ Fn×m

2

, we have QP(f ) − I(f ) = 1 2m ∑

x∈P

∑

k∈Ns

ˆ f (k)walk(x) − ˆ f (0)

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 36 / 61

slide-107
SLIDE 107

Decomposition of integration error

For a digital net P with C1, . . . , Cs ∈ Fn×m

2

, we have QP(f ) − I(f ) = 1 2m ∑

x∈P

∑

k∈Ns

ˆ f (k)walk(x) − ˆ f (0) = ∑

k∈Ns

ˆ f (k) 1 2m ∑

x∈P

walk(x) − ˆ f (0)

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 36 / 61

slide-108
SLIDE 108

Decomposition of integration error

For a digital net P with C1, . . . , Cs ∈ Fn×m

2

, we have QP(f ) − I(f ) = 1 2m ∑

x∈P

∑

k∈Ns

ˆ f (k)walk(x) − ˆ f (0) = ∑

k∈Ns

ˆ f (k) 1 2m ∑

x∈P

walk(x) − ˆ f (0) = ∑

k∈P⊥\{0}

ˆ f (k)

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 36 / 61

slide-109
SLIDE 109

Decomposition of integration error

For a digital net P with C1, . . . , Cs ∈ Fn×m

2

, we have QP(f ) − I(f ) = 1 2m ∑

x∈P

∑

k∈Ns

ˆ f (k)walk(x) − ˆ f (0) = ∑

k∈Ns

ˆ f (k) 1 2m ∑

x∈P

walk(x) − ˆ f (0) = ∑

k∈P⊥\{0}

ˆ f (k) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + ∑

k∈P⊥\{0} 2n|k

ˆ f (k)

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 36 / 61

slide-110
SLIDE 110

Decomposition of integration error

For a digital net P with C1, . . . , Cs ∈ Fn×m

2

, we have QP(f ) − I(f ) = 1 2m ∑

x∈P

∑

k∈Ns

ˆ f (k)walk(x) − ˆ f (0) = ∑

k∈Ns

ˆ f (k) 1 2m ∑

x∈P

walk(x) − ˆ f (0) = ∑

k∈P⊥\{0}

ˆ f (k) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + ∑

k∈P⊥\{0} 2n|k

ˆ f (k) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + ∑

k∈Ns

0\{0}

2n|k

ˆ f (k).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 36 / 61

slide-111
SLIDE 111

What is the second term?

QP(f ) − I(f ) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + ∑

k∈Ns

0\{0}

2n|k

ˆ f (k) The second term is...

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 37 / 61

slide-112
SLIDE 112

What is the second term?

QP(f ) − I(f ) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + ∑

k∈Ns

0\{0}

2n|k

ˆ f (k) The second term is... ∑

k∈Ns

0\{0}

2n|k

ˆ f (k) = ∑

k∈Ns

ˆ f (k) 1 2ns

2n−1

∑

i1,...,is=0

walk ( i1 2n , . . . , is 2n ) − ˆ f (0) = 1 2ns

2n−1

∑

i1,...,is=0

∑

k∈Ns

ˆ f (k)walk ( i1 2n , . . . , is 2n ) − I(f ) = 1 2ns

2n−1

∑

i1,...,is=0

f ( i1 2n , . . . , is 2n ) − I(f ) the integration error for a regular grid!

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 37 / 61

slide-113
SLIDE 113

Digital net as a subset of regular grid

QP(f ) − I(f ) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) +   1 2ns

2n−1

∑

i1,...,is=0

f ( i1 2n , . . . , is 2n ) − I(f )  

Figure: Hammersley point set for m = n = 6 with regular grid

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 38 / 61

slide-114
SLIDE 114

Euler-Maclaurin

If f has partial mixed derivatives up to order α in each variable, 1 2ns

2n−1

∑

i1,...,is=0

f ( i1 2n , . . . , is 2n ) − I(f ) = c1(f ) 2n + · · · + cα−1(f ) 2(α−1)n + O(2−αn), where cτ(f ) = ∑

τ1,...,τs≥0 τ1+···+τs=τ

I(f (τ1,...,τs))

s

∏

j=1 τj>0

bτj with bi denoting the i-th Bernoulli number divided by i!. For a detailed proof of this claim, see Dick, G. & Yoshiki (2019).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 39 / 61

slide-115
SLIDE 115

Euler-Maclaurin

Using this result, we have QP(f ) − I(f ) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + c1(f ) 2n + · · · + cα−1(f ) 2(α−1)n + O(2−αn).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 40 / 61

slide-116
SLIDE 116

Euler-Maclaurin

Using this result, we have QP(f ) − I(f ) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + c1(f ) 2n + · · · + cα−1(f ) 2(α−1)n + O(2−αn). Consider high precision, say n = αm. Noting that N = 2m, the first term becomes dominant since QP(f ) − I(f ) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + c1(f ) Nα + · · · + cα−1(f ) Nα(α−1) + O(N−α2). Remark: For functions with dominating mixed smoothness α, the integration error cannot be better than O(N−α) in general.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 40 / 61

slide-117
SLIDE 117

Richardson extrapolation

Consider the case n = m. Let us simplify the situation for a moment and consider the following expansions for successive values of m:

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 41 / 61

slide-118
SLIDE 118

Richardson extrapolation

Consider the case n = m. Let us simplify the situation for a moment and consider the following expansions for successive values of m: I (1)

m−α+1 := c0 +

c1 2m−α+1 + c2 22(m−α+1) + · · · + cα−1 2(α−1)(m−α+1) . . . I (1)

m−1 := c0 +

c1 2m−1 + c2 22(m−1) + · · · + cα−1 2(α−1)(m−1) I (1)

m

:= c0 + c1 2m + c2 22m + · · · + cα−1 2(α−1)m Problem: Estimate c0 from values of I (0)

m−α+1, . . . , I (0) m

as precisely as possible (without knowing c1, . . . , cα−1).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 41 / 61

slide-119
SLIDE 119

Richardson extrapolation

Do the following recursion: For τ = 1, . . . , α − 1, compute I (τ+1)

n

:= 2τI (τ)

n

− I (τ)

n−1

2τ − 1 for n = m − α + τ + 1, . . . , m Schematically this recursion goes like

I (1)

m−α+1

I (1)

m−α+2

· · · I (1)

m−1

I (1)

m

↘ ↓ ↘ · · · ↓ ↘ ↓ I (2)

m−α+2

· · · I (2)

m−1

I (2)

m

↘ · · · ↓ ↘ ↓ ... . . . I (α−1)

m−1

I (α−1)

m

↘ ↓ I (α)

m

We will end up with I (α)

m

= c0. In reality, we have other terms...

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 42 / 61

slide-120
SLIDE 120

Use a set of digital nets

For m − α < n ≤ m, let Pn be a digital net with C1, . . . , Cs ∈ Fn×n

2

. Here we do NOT require Pm−α+1 ⊂ Pm−α+2 ⊂ · · · ⊂ Pm. Richardson extrapolation can push out the intermediate terms of the expansion QPn(f ) − I(f ) = ∑

k∈P⊥

n \{0}

2n∤k

ˆ f (k) + c1(f ) 2n + · · · + cα−1(f ) 2(α−1)n + O(2−αn) without knowing c1(f ), . . . , cα−1(f ).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 43 / 61

slide-121
SLIDE 121

Use a set of digital nets

This way we arrive at an extrapolated QMC rule which satisfies

m

∑

n=m−α+1

wnQPn(f ) − I(f ) =

m

∑

n=m−α+1

wn ∑

k∈P⊥

n \{0}

2n∤k

ˆ f (k) + O(2−αm), where wn’s depend only on α,

m

∑

n=m−α+1

wn = 1 and

m

∑

n=m−α+1

|wn| ≤

α−1

∏

i=1

2i + 1 2i − 1. NOTE: The total number of points is 2m−α+1 + · · · + 2m < 2m+1.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 44 / 61

slide-122
SLIDE 122

What we want

m

∑

n=m−α+1

wnQPn(f ) − I(f ) =

m

∑

n=m−α+1

wn ∑

k∈P⊥

n \{0}

2n∤k

ˆ f (k) + O(2−αm), If we have digital nets with C1, . . . , Cs ∈ Fn×n

2

which satisfy ∑

k∈P⊥

n \{0}

2n∤k

ˆ f (k) = O(2−(α−ε)n) with arbitrarily small ε > 0, extrapolated QMC rule achieves the almost optimal rate of convergence for smooth functions.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 45 / 61

slide-123
SLIDE 123

Outline

QMC and digital nets/sequences Classical QMC Higher order QMC Richardson extrapolation and QMC Application 1: Truncation of higher order nets and sequences Application 2: Extrapolated polynomial lattice rules

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 46 / 61

slide-124
SLIDE 124

Truncation operator

For x = (0.ξ1ξ2 . . .)2 ∈ [0, 1), define the map trm(x) := (0.ξ1ξ2 . . . ξm00 . . .)2, For a point set P ⊂ [0, 1)s, we write trm(P) = {trm(x) | x ∈ P}, where trm is applied componentwise. If P is a digital net with C1, . . . , Cs ∈ Fn×m

2

for n ≥ m, trm(P) is a digital net with upper m × m submatrices of C1, . . . , Cs.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 47 / 61

slide-125
SLIDE 125

What G. (2019+) proves

Theorem

Let P be an order α digital (t, m, s)-net. For functions with dominating mixed smoothness α ≥ 2, we have ∑

k∈(trm(P))⊥\{0} 2m∤k

ˆ f (k) = O ((log N)αs Nα ) , with N = 2m. Compare: ∑

k∈P⊥\{0}

ˆ f (k) = O ((log N)αs Nα ) but ∑

k∈(trm(P))⊥\{0}

ˆ f (k) = O ( 1 N )

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 48 / 61

slide-126
SLIDE 126

Truncation of higher order digital nets

Figure: Interlaced Sobol’ point sets: original (left) and after truncation (right)

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 49 / 61

slide-127
SLIDE 127

Numerical results: original higher order QMC

Univariate test function: f (x) = x3 (log x + 0.25) . Interlaced Sobol’ sequences: If αm > 52, tr52 is applied.

5 10 15 20 25 log

2N

  • 70
  • 60
  • 50
  • 40
  • 30
  • 20
  • 10

log

2(absolute integration error) Order 2 Sobol' Order 3 Sobol' N-2 N-3

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 50 / 61

slide-128
SLIDE 128

Numerical results: extrapolation

QP(f ) − I(f ) = ∑

k∈P⊥\{0} 2n∤k

ˆ f (k) + c1(f ) 2n + · · · + cα−1(f ) 2(α−1)n + O(2−αn). Extrapolation applied to truncated interlaced Sobol’ sequences: order 2 (left) and order 3 (right)

5 10 15 20 25 log

2N

  • 50
  • 45
  • 40
  • 35
  • 30
  • 25
  • 20
  • 15
  • 10
  • 5

log

2(absolute integration error) Im (1) Im (2) N-1 N-2

5 10 15 20 25 log

2N

  • 70
  • 60
  • 50
  • 40
  • 30
  • 20
  • 10

log

2(absolute integration error) Im (1) Im (2) Im (3) N-1 N-2 N-3

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 51 / 61

slide-129
SLIDE 129

Numerical results: comparison for large s

Test function: f (x) =

s

∏

j=1

[ 1 + γj ( xc

j −

1 c + 1 )] , with s = 100, c = 1.3 and γj = j−2.

5 10 15 20 log

2N

  • 40
  • 35
  • 30
  • 25
  • 20
  • 15
  • 10
  • 5

log

2(absolute integration error) Im (2) Jm (2) Order 2 Sobol' N-2

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 52 / 61

slide-130
SLIDE 130

Outline

QMC and digital nets/sequences Classical QMC Higher order QMC Richardson extrapolation and QMC Application 1: Truncation of higher order nets and sequences Application 2: Extrapolated polynomial lattice rules

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 53 / 61

slide-131
SLIDE 131

What Dick, G. & Yoshiki (2019) prove

Theorem

Let p ∈ F2[x] be irreducible with deg(p) = m. Component-by-component algorithm based on a computable criterion can find good q ∈ (F2[x])s which satisfies ∑

k∈P(p,q)⊥\{0} 2m∤k

ˆ f (k) = O(2−(α−ε)m) with arbitrarily small ε > 0, for functions with dominating mixed smoothness α. Here “a computable criterion” has been provided by Baldeaux, Dick, Leobacher, Nuyens & Pillichshammer (2012), which stems from the decay of Walsh coefficients for smooth functions.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 54 / 61

slide-132
SLIDE 132

Construction cost

Since we do NOT require Pm−α+1 ⊂ Pm−α+2 ⊂ · · · ⊂ Pm, polynomial lattice point set of each size can be constructed independently. Construction cost of each is of O((s + α)n2n). In total, the cost will be of order

m

∑

n=m−α+1

(s + α)n2n ≤ (s + α)N log2 N, where N = 2m−α+1 + · · · + 2m. This is better than interlaced polynomial lattice point sets for which we need O(αsN log N).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 55 / 61

slide-133
SLIDE 133

Fast QMC matrix-vector product

In some applications, we want to estimate the integral of the form ∫

[0,1)s f (Ax) dx.

In terms of computational cost, computing Ax can be dominant. Even if A does not have any special structure, this computation can be made fast by using (polynomial) lattice point sets (Dick, Kuo, Le Gia, Schwab, 2015). This strategy does not work for interlaced polynomial lattice rules, whereas it does work for extrapolated polynomial lattice rules.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 56 / 61

slide-134
SLIDE 134

Numerical results (α = 2)

Test function: f (x) =

s

∏

j=1

[ 1 + γj 1 + γjx2

j

] , with s = 100 and γj = j−2. Extrapolated rule with α = 2 (red), interlaced rule with α = 2 (green)

4 6 8 10 12 14 16 log

2N

10-12 10-10 10-8 10-6 10-4 10-2 100

Absolute error

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 57 / 61

slide-135
SLIDE 135

Numerical results (α = 3)

Test function: f (x) =

s

∏

j=1

[ 1 + γj 1 + γjx2

j

] , with s = 100 and γj = j−2. Extrapolated rule with α = 3 (red), interlaced rule with α = 3 (green)

4 6 8 10 12 14 16 log

2N

10-12 10-10 10-8 10-6 10-4 10-2 100

Absolute error

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 58 / 61

slide-136
SLIDE 136

Conclusions

As far as I know, this is the first attempt to apply an extrapolation technique to QMC. Roughly speaking, Richardson extrapolation is shown useful for reducing necessary precision significantly. Two results are obtained so far

1

truncation of higher order digital nets/sequences: high-order convergence can be recovered by applying Richardson extrapolation.

2

introduction of extrapolated polynomial lattice rules: without relying on interlacing, we can construct high-order convergent QMC-based rules.

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 59 / 61

slide-137
SLIDE 137

Conclusions

As far as I know, this is the first attempt to apply an extrapolation technique to QMC. Roughly speaking, Richardson extrapolation is shown useful for reducing necessary precision significantly. Two results are obtained so far

1

truncation of higher order digital nets/sequences: high-order convergence can be recovered by applying Richardson extrapolation.

2

introduction of extrapolated polynomial lattice rules: without relying on interlacing, we can construct high-order convergent QMC-based rules.

Do we have another good application of extrapolation? What about lattice rules? Can we randomize extrapolated QMC rule? Partially yes, if the same randomly chosen digital shift is applied. What about scrambling?

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 59 / 61

slide-138
SLIDE 138

If you are interested, please refer to

  • J. Dick, TG, T. Yoshiki: Richardson extrapolation of polynomial

lattice rules, SIAM Journal on Numerical Analysis, 57(1):44-69, 2019. TG: Richardson extrapolation allows truncation of higher order digital nets and sequences, IMA Journal of Numerical Analysis, online (arXiv:1805.09490).

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 60 / 61

slide-139
SLIDE 139

Thank you for your attention!

Takashi Goda (U. Tokyo) Extrapolation and HOQMC MCM 2019 61 / 61