Unidirectional principle for non-local operators

The unidirectional principle for differential operators [BZ96, Bun92b, Zei11] is a sophisticated algorithm that realizes the matrix-vector multiplication for generalized sparse grid problems in a computational complexity that is strictly linear in the number of degrees of freedom. Now, we present a generalization from [GH13c], where it is no longer necessary to specifically tailor the TopDown/BottomUp algorithms to the operator in use. Moreover, it can be employed with non-local operators like integro-differential operators. A related abstraction of the unidirectional principle for multilevel discretizations has been presented in [Zei11].

With the same assumptions as in the previous subsection, we want to computey_I =A_Ix_I,

Algorithm 1Propagation phase of single space matrix-vector multiplication Input: multi index k, index set I and subspace vector xk

Output: z_l,k= N

p:lp<kpI^(p)_l

p,kpA^(p)_k

p,kp⊗N

q:lq≥kqI^(q)_l

q,kq

x_k,l∈ I

1: Mark eachlinI as “not visited”

2: z_k,k←x_k .Start with subspace vector k

3: Markkas “visited”

4: for all l∈ I do

5: if l is marked as “not visited” then

6: computeSubspaceVec(l)

7: end if

8: end for

9: returnz_l,k,l∈ I

10: functioncomputeSubspaceVec(l)

11: if l is marked as “visited” then

12: return

13: end if

14: forp= 1, . . . , d do .Since l6=k, there is at least one pwithl_p 6=k_p

15: if lp > kp then

16: l⁰ ←l−e_p .Decrease p-th entry by 1

17: computeSubspaceVec(l’) . z_l⁰_,k is now set

18: z_l,k← I^(p)_l

p,lp−1⊗N

q6=pI^(q)_l

q,lq

z_l⁰_,k . Prolongate from lower subspace

19: Markl as “visited”

20: return

21: else if l_p =k_p−1 then

22: l⁰ ←l+ep .Increase p-th entry by 1

23: computeSubspaceVec(l’) . z_l⁰_,k is now set

24: z_l,k← (I^(p)_l

p,kpA^(p)_k

p,kp)⊗N

q6=pI^(q)_l

q,lq

z_l⁰_,k . Apply matrixA^(p)_kp,kp and restrict result 25: Marklas “visited”

26: return

27: else if lp< kp−1then

28: l⁰←l+ep .Increasep-th entry by1

29: computeSubspaceVec(l’) .zl⁰,k is now set 30: zl,k← I^(p)_lp,lp+1⊗N

q6=pI^(q)_lq,lq

zl⁰,k .Restrict result from higher subspace 31: Marklas “visited”

32: return

33: end if 34: end for 35: end function

4.3 Operator application 63 Algorithm 2Application phase of single space matrix-vector multiplication

Input: multi indexk, index set I andz_l,k= N

p:lp<kpI^(p)_l

p,kpA^(p)_k

p,kp⊗N

q:lq≥kqI^(q)_l

q,kq

x_k,l∈ I Output: w_l,k=

N_d

p=1A^(p)_l

p,kp

x_k,l∈ I

1: for all linI do

2: w_l,k←

N

p:lp<kpI^(p)_l

p,lp⊗N

q:lq≥kqA^(q)_l

q,lq

z_l,k

3: end for

i.e.,

y_l =X

k∈I

A⁽¹⁾_l

1,k1⊗ · · · ⊗A^(d)_l

d,kp

x_k (4.37)

for alll∈ I.

First, we start with the one-dimensional case ofI=F_J¹ and split the sum (4.37) in two parts y_l =

J

X

k=0

A_l,kx_k

=

l

X

k=0

A_l,kx_k+

k

X

k=l+1

A_l,kx_k . (4.38)

Both sums in (4.38) need to be treated separately. As already indicated, the representa-tion (4.18) of prolongarepresenta-tions and restricrepresenta-tions is a very efficient tool to transport intermediate results. The computation of the first sum in (4.38) can be done using the TopDown algorithm.

TopDown algorithm

The TopDown algorithm uses a simple recursive relation to compute y¯_l:=

l

X

k=0

A_l,kx_k (4.39)

for all0≤l≤J. First, we define vectors z_l:=

l

X

k=0

I_l,kx_k, for l= 0, . . . , J that satisfy the recursive relationship

z_l=

l

X

k=0

I_l,kx_k=

l−1

X

k=0

I_l,l−1I_l−1,kx_k+x_l =I_l,l−1z_l−1+x_l

for l≥1. Now, (4.39) can be expressed by

¯ y_l=

l

X

k=0

A_l,lI_l,kx_k=A_l,l

l

X

k=0

I_l,kx_k=A_l,lz_l . (4.40) As long as the prolongationsI_l,l−1 and the application of the matricesA_l,l work in linear time, the TopDown algorithm is of optimal order. The whole procedure is described in pseudocode in Algorithm 3.

Algorithm 3TopDown algorithm

Input: MatricesA_l,l and x_l for l= 0, . . . , J Output: y¯_l=Pl

k=0A_l,kx_k, l= 0, . . . , J

1: z₀ ←x₀

2: y¯₀←A_0,0z₀

3: for all l= 1, . . . , J do

4: z_l ←I_l,l−1z_l−1+x_l

5: y¯_l←A_l,lz_l

6: end for

BottomUp algorithm

The BottomUp algorithm computes

y_l:=

J

X

k=l+1

A_l,kx_k (4.41)

for alll= 0, . . . , J−1, and we definey_J :=0. The recursive relationship

yl =

J

X

k=l+2

A_l,kx_k+A_l,l+1x_l+1

= I_l,l+1 X^J

k=(l+1)+1

A_l+1,kx_k+A_l+1,l+1x_l+1

= I_l,l+1

yl+1+A_l+1,l+1x_l+1

. (4.42)

holds for l=J−1, . . . ,0. Clearly, ally

l can be precalculated in linear time provided that the restrictionsI_l,l+1 and the application of the matricesA_l,l work in linear time. The pseudocode is given in Algorithm 4.

4.3 Operator application 65 Algorithm 4BottomUp algorithm

Input: MatricesA_l,l and x_l for l= 0, . . . , J Output: y_l=P_J

k=l+1A_l,kx_k, l= 0, . . . , J

1: y

J ←0

2: for all l=J−1, . . . ,0do

3: y

l←I_l,l+1

yl+1+A_l+1,l+1x_l+1

4: end for

Multi-dimensional case

The multi-dimensional case can be reduced to the recursive application of the one-dimensional algorithms TopDown and BottomUp. To this end, we split and rearrange the sum

y_l =X

k∈I

A⁽¹⁾_l

1,k1⊗A⁽²⁾_l

2,k2 ⊗ · · · ⊗A^(d)_l

d,kd

x_k ∀l∈ I . So, for all l∈ I, we need to compute

yl = X

k1≤l1with (k1,l2,...,ld)∈I

(A⁽¹⁾_l

1,k1 ⊗I⁽²⁾_l

2,l2⊗ · · · ⊗I^(d)_l

d,ld)· (4.43)

X

(k2,...,kd)with (k1,k2,...,kd)∈I

I⁽¹⁾_k

1,k1⊗A⁽²⁾_l

2,k2 ⊗ · · · ⊗A^(d)_l

d,kd

x_{(k1,k2,...,kd)} (4.44)

+ X

(k2,...,kd)with (l1,k2,...,kd)∈I

I⁽¹⁾_l

1,l1 ⊗A⁽²⁾_l

2,k2⊗ · · · ⊗A^(d)_l

d,kd

· (4.45)

X

k1>l1with (k1,k2,...,kd)∈I

(A⁽¹⁾_l

1,k1 ⊗I⁽²⁾_k

2,k2 ⊗ · · · ⊗I^(d)_k

d,kd)x_{(k1,k2,...,kd)}. (4.46) Now, (4.43) resembles the application of the one-dimensional TopDown algorithm, and (4.46) re-sembles the application of the one-dimensional BottomUp algorithm. The sums (4.44) and (4.45) are the result of a recursive application of the multi-dimensional algorithm with the first dimen-sion left unchanged. The monotonicity condition (4.8) ensures that intermediate results can be stored with the same complexity asx_I and the final result y_I. Specifically, we know that

(l₁, . . . , l_d)∈ I, k₁ ≤l₁ ⇒(k₁, l₂, . . . , l_d)∈ I in (4.43) and

(k₁, . . . , k_d)∈ I, k₁ > l₁⇒(l₁, k₂, . . . , k_d)∈ I

in (4.45), which means that intermediate results can be represented as generalized sparse grid functions with the same index set I. See Algorithm 5 for a description in pseudocode, where we use the shorthand notation

k⁰∪ {k_p}= (k₁, . . . , k_d)

for k⁰ = (k1, . . . , kp−1, kp+1, . . . , kd)∈N^d−1. Algorithm 5Unidirectional principle

Input: Index setI, Matrices A^(p)_l,k,l,k∈ I and p= 1, . . . , d, vector x_I Output: y_I = (yl)l∈I withyl =P

k∈I A⁽¹⁾_l

1,k1⊗ · · · ⊗A^(d)_l

d,kp

xk for alll∈ I

1: y_I ←UniDir(0,x_I)

2: functionUniDir(p,x_I) . pis active dimension, x_I the input vector

3: if p=dthen

4: return . Nothing to do

5: end if

6: z_I ←x_I . Copy input vector

7: if p6= 0 then

8: for all k⁰ ∈N^d−1 for which ∃k_p :k⁰∪ {k_p} ∈ I do

9: BottomUp((x_k)_k:=k⁰_∪{kp}∈I,(A^(p)_k

p,kp)_kp:k⁰_∪{kp}∈I). Apply (4.46) tox_I in place

10: end for

11: end if

12: UniDir(p+ 1,x_I) . Recursive call (4.45)

13: UniDir(p+ 1,z_I) . Recursive call (4.44)

14: if p6= 0 then

15: for all k⁰ ∈N^d−1 for which ∃k_p :k⁰∪ {k_p} ∈ I do

16: TopDown((z_k)_k:=k⁰_∪{kp}∈I,(A^(p)_k

p,kp)_kp:k⁰_∪{kp}∈I) .Apply (4.43) toz_I in place

17: end for

18: end if

19: x_I ←x_I+z_I

20: end function

Note here that the cost complexity of the algorithm is only linear with respect to the degrees of freedom if the cost complexity of the univariate operator applications is linear. However, two recursive calls per dimension lead to 2^d calls of the function. This exponential dependence of the computational complexity on the dimension is undesirable, but can be avoided in special cases, see [Feu05]. Note that the parallelization of the unidirectional principle is far from trivial due to the involved recursive nature of the algorithm, see [HSB12].

5 Preconditioning of high-dimensional elliptic equations

In this chapter, we develop additive (sometimes coined “parallel”) Schwarz preconditioners for generalized sparse grid discretizations of symmetricH-elliptic variational problems

a(u, v) =F(v) ∀v∈H , (5.1)

whereHis thek·ka:=a(·,·)¹²-closure ofspan{Vl:l∈N^d}andF is a bounded linear functional on H. The main focus of this chapter is on generic H^t- orH_mix^t,l -elliptic problems. Ultimately, we will apply our preconditioner to variational problems (3.16) that arise in the course of solving the BKE, even though they can be asymmetric. Our preconditioners still work well in these cases. Note that by our choice of subspaces in the additive Schwarz setting, we obtain methods that are in fact preconditioners for the system of linear equations (4.22), which is the generalized sparse grid version of (3.19).

For isotropic full grid discretizations, an optimal iteration count which is independent of the number of degrees of freedom is typically achieved by multiplicative multigrid methods [Yse93, BL11, Hac85, Gri94b], the additive BPX preconditioner [BPX90, Osw92, Osw94] or wavelet-based methods. In the sparse grid case, the condition number is more difficult to reduce than in the regular full grid case. For example, already for a straightforward regular sparse grid discretization, cf. [GO94], a simple diagonal scaling similar to the case of the BPX-preconditioner does not result in asymptotically bounded condition numbers in dimensions d≥3. Here, more complicated basis functions like prewavelets offer a solution [GO95b].

The main idea is to apply a subspace correction method to (5.1), where the subspaces solvers are based on auxiliary bilinear forms on the anisotropic full grid spaces that compose the sparse grid. Their relative scaling is at our disposal and amounts to a diagonal scaling of the operator matrix of the discretized system (4.22). We heavily rely on a norm equivalence based on orthogonal complement spaces. First we show that, for a certain class of second-order PDEs, the norm equivalence constants are independent of the space dimension and the diffusion coefficients. Then, based on the norm equivalence, we infer quasi-optimal scaling factors for our full grid spaces by a Linear Program (LP). This approach closely follows the lines of [GHO15].

We prove that O(J^d−2) is a lower bound for the condition number for any positive diagonal scaling.

This motivates the use of partially negative scaling factors, and we present an algebraic transformation that results in optimal condition numbers. In fact, we even observe falling condition numbers with rising dimension for the Laplace operator discretized by a regular sparse grid. We also present a method closely related to a prewavelet approach which results in exactly the same condition numbers as the algebraic transformation, but only needs positive scalings and produces symmetrizable matrices, see also [GH14b].

67

If a norm equivalence is not available, we can still employ a non-linear variable precondi-tioner [JN99] that has been referred to as OptiCom in the case of sparse grids in data min-ing [Heg03, Gar06, HGC07]. We describe an efficient implementation based on the smin-ingle space matrix-vector multiplication presented in Chapter 4. This reduces the typically quadratic costs of the OptiCom to log-linear with respect to the degrees of freedom if the associated bilinear form a(·,·) is given as a sum of tensor products. All preconditioners allow a CG version, and with the exception of the variable preconditioner, possess a cost complexity that is only linear with respect to the degrees of freedom.

This chapter is organized as follows: In Section 5.1, we describe a well-known equivalence of norms for H^t-elliptic problems and prove that its norm equivalence constants are independent of the diffusion coefficients and the space dimension. In Section 5.2, we build the connection between the norm equivalence and our generalized sparse grid space via the theory of subspace splittings. In Section 5.3, we derive a Linear Program to find close-to-optimal positive scalings but prove the limits of this approach forH^t-elliptic problems. In Section 5.4, we admit negative values for the scaling of our generating system and thus obtain optimal condition numbers.

The same condition numbers are realized by a block-diagonal preconditioner in Section 5.5. In Section 5.6, the non-linear variable preconditioner OptiCom is described. Section 5.7 contains experiments up to dimensiond= 10that support our theoretical findings.

5.1 Norm equivalences based on orthogonal subspaces

Introducing the orthogonal complement spaces W_l^(p) = V_l^(p) _V(p) V_l−1^(p) for l ≥ 1, and setting W₀^(p) := V₀^(p), we can writeV as the (·,·)V-orthogonal (from now on V-orthogonal) sumV = L

l∈N^dW_l of the subspaces W_l=W_l⁽¹⁾

1 ⊗. . .⊗W_l^(d)

d , l= (l₁, . . . , l_d)∈N^d. Similarly, any finite-dimensional space V_l =V_l⁽¹⁾

1 ⊗. . .⊗V_l^(d)

d ⊂V,l∈N^d, is theV-orthogonal sum V_l = L

k≤lW_k, where the inequality k ≤ l is meant componentwise. Furthermore, it is easy to verify that

V_I =M

k∈I

W_k (5.2)

can be written as V-orthogonal sum of the subspacesW_k withk∈ I.

Based on (5.2), we can estimate the redundancy we create by using a generating system instead of a basis. Let us assume that for our spaces V_l^(p), l∈N, p= 1, . . . , d,

c·n_l−1 ≤n_l ⇔ n_l−1 ≤ ¹_c ·n_l holds forl≥1, e.g.,c= 2for dyadically refined linear splines. Then,

dimW_l=n_l−n_l−1≥(1−¹_c)n_l

5.1 Norm equivalences based on orthogonal subspaces 69 and

dimV_I =X

k∈I

dimWk≥(1−¹_c)^dX

k∈I

nk ⇒N_I ≤(1−¹_c)^−ddimV_I , (5.3) which means we can limit the number of degrees of freedom in our generating systemN_I by the dimension ofV_I times a factor(1−¹_c)^−d. It is not hard to see that this factor is relatively sharp, and thus, for c = 2, we have a factor of about 2^d more degrees of freedom in our generating system than needed in the sparse grid space. This result coincides with the somewhat simpler calculation (4.17).

The basis for many further considerations in this thesis is the assumption of a set of fixed positive weightsβk,k∈N^d,and a norm equivalence

λ_min X

k∈N^d

β_kkw_kk²V ≤ kuk²a≤λ_max X

k∈N^d

β_kkw_kk²V , (5.4) where0 < λ_min ≤λ_max<∞ and the w_k ∈W_k,k∈N^d,denote the components of the unique V-orthogonal decomposition of u∈ H, i.e., u =P

k∈N^dw_k. The estimate (5.4) is assumed to be sharp, which means that the norm equivalence constants λ_min andλ_maxare given by

λ_min := inf

06=u∈H

kuk²a

P

k∈N^dβ_kkw_kk²V

and λ_max:= sup

06=u∈H

kuk²a

P

k∈N^dβ_kkw_kk²V

.

From now on, we use the symbol ' to indicate such an equivalence, and call κ = λ_max/λ_min condition number of the splitting.

The described setup is motivated by the discretization ofH^t-elliptic problems by generalized sparse grid spaces V_I over thed-dimensional unit cube. The restriction ont is|t|< r+ 3/2 if we useC^r spline spaces of fixed degreem≥r+ 1 over dyadic partitions of step size2^−l as the building blocks V_l^(p). Then, the equivalence of norms (5.4) has the form

kuk²_H^t ' X

k∈N^d

2^2t|k|∞kwkk²_L2(Ω^d₎, |k|∞= max

i=1,...,dk_i , (5.5)

where the wk ∈ Wk,k ∈ N^d, denote the L₂-orthogonal components of the function u ∈ H^t, see [Osw94] for this kind of results. A more general result for H_mix^t,l -elliptic problems is stated in [GK09].

We now single out the case of linear spline spaces for H¹-elliptic problems and take a look at the norm equivalence constants for a special set of weights (β_k)_k∈Nd. Due to Jackson- and Bernstein-inequalities [Dah96, Osw94] for dyadically refined linear spline spaces (V_l^(p))^∞_l=1, p= 1, . . . , d, we have a norm equivalence

∂u

∂x,∂u

∂x

L2(Ω)

'X

l∈N

β_lkw_lk²L2(Ω) (5.6) foru∈H₀¹(Ω)withβ_l= 2^2l,u=P

l∈Nw_l, wherew_l∈W_l, l∈N, and0< λ⁽¹⁾_min≤λ⁽¹⁾_max<∞. Of course, (5.6) holds also foru∈V_J^(p), p= 1, . . . , d,with norm equivalence constantsλ^(1),J_min ≥λ⁽¹⁾_min and λ^(1),J_max ≤ λ⁽¹⁾_max for all J ∈ N. The next theorem shows that this norm equivalence carries

over to the d-dimensional case with exactly the same norm equivalence constants. This is a generalization of the proof given in [GH14b].

Theorem 5.1. For u∈H₀¹(Ω^d), αp>0, p= 1, . . . , d and a(u, v) =

d

X

p=1

αp

∂u

∂x_p, ∂v

∂x_p

L2(Ω^d)

, (5.7)

it holds that a(u, u)' X

l∈N^d

X^d

p=1

αp2^2lp

kwlk²_L2(Ω^d₎ for u= X

l∈N^d

wl with wl∈Wl,l∈N^d, (5.8)

where the constantsλ^(d)_minandλ^(d)maxassociated with(5.8)are the same as in (5.6)i.e. λ^(d)_min=λ⁽¹⁾_min and λ^(d)max = λ⁽¹⁾max. For functions u ∈ V_Fd

J = VJ ⊂ H₀¹(Ω^d), the corresponding equivalence constantsλ^(d),Jmax and λ^(d),J_min of (5.8) are equal to λ^(1),Jmax and λ^(1),J_min , respectively, forJ ∈N.

Proof. We start with the case u ∈ V_F^d

J, J ∈ N. In Chapter 4, the functions (φ_J,i)ⁿ_i=1^J formed a basis for the spaces V_J^(p), p= 1, . . . , d. Of course, there also exists an L2-orthonormal basis (ψ_J,i)ⁿ_i=1^J of V_J. Furthermore, we need the orthogonal decomposition(ω_l,i)^J_l=1 of ψ_J,i∈ V_J for all i= 1, . . . , n_J withω_l,i∈W_l, l= 1, . . . , J and

ψ_J,i=

J

X

l=1

ω_l,i. Next, analogously to (4.5), we define

ψ_J,i(x) =ψ_J,i1(x₁)· · ·ψ_J,id(x_d) and ω_l,i(x) =ω_l1,i1(x₁)· · ·ω_ld,id(x_d)

for x= (x₁, . . . , x_d),i ∈χJ and l∈ F_J^d. This opens a direct way to find orthogonal decompo-sitions of functionsu=P

i∈χJz_iψ_J,i∈V_Fd J by

u= X

i∈χJ

z_i X

l∈FJ^d

ω_l,i= X

l∈FJ^d

X

i∈χJ

z_iω_l,i= X

l∈FJ^d

w_l

with

wl= X

i∈χJ

ziω_l,i∈Wl (5.9)

for alll∈ FJ^d.

Now, we show that the norm equivalence (5.8) holds for any u ∈ V_Fd

J with the constants

5.1 Norm equivalences based on orthogonal subspaces 71 λ^(1),J_max andλ^(1),J_min from (5.6). It holds that

a(u, u) =

d

X

p=1

α_p ∂

∂x_p X

i∈χJ

z_iψ_J,i, ∂

∂x_p X

j∈χJ

z_jψ_J,j

L2(Ω^d) (5.10)

=

d

X

p=1

α_pX

i∈χJ

X

j∈χJ

∂

∂xp

z_iψ_J,ip, ∂

∂xp

z_jψ_J,jp

L2(Ω) d

Y

q=1 q6=p

(ψ_J,iq, ψ_J,jq)_L2(Ω) (5.11)

=

d

X

p=1

X

i⁰=i\{ip} i∈χJ

∂

∂xp nJ

X

ip=1

z_i⁰_∪{ip}ψ_J,ip, ∂

∂xp nJ

X

jp=1

z_i⁰_∪{jp}ψ_J,jp

L2(Ω). (5.12) We obtain (5.11) by repeated applications of the distributive law and by using the product structure of the L₂-scalar product. Then, the orthonormal basis property of the (ψ_J,i)ⁿ_i=1^J cancels all terms fori_q6=j_q, q6=p and we get (5.12). Here we use the notation

i \ {i_p}= (i₁, . . . , i_p−1, i_p+1, . . . , i_d) and i⁰ ∪ {i_p}= (i₁, . . . , i_p−1, i_pi_p+1, . . . , i_d).

We can apply the one-dimensional norm equivalence (5.6) to (5.12) and obtain the upper bound

· · · ≤

d

X

p=1

α_pλ^(1),J_max X

i⁰=i\{ip} i∈χJ

J

X

lp=1

2^2lpX^nJ

ip=1

z_i⁰_∪{ip}ω_lp,ip,

nJ

X

jp=1

z_i⁰_∪{jp}ω_lp,jp

L2(Ω) (5.13)

=λ^(1),J_max

d

X

p=1

α_p X

i∈χJ

X

j∈χJ

J

X

lp=1

2^2lp(z_iω_lp,ip, z_jω_lp,jp)_L2(Ω)·

d

Y

q=1 q6=p

(ψ_J,iq, ψ_J,jq)_L2(Ω) (5.14)

=λ^(1),J_max

d

X

p=1

α_p X

i∈χJ

X

j∈χJ

X

l∈F_J^d

2^2lp(z_iω_lp,ip, z_jω_lp,jp)_L2(Ω)·

d

Y

q=1 q6=p

(ω_lq,iq, ω_lq,jq)_L2(Ω) (5.15)

=λ^(1),J_max X

l∈F_J^d

X^d

p=1

α_p2^2lp X

i∈χJ

z_iω_l,i,X

j∈χJ

z_jω_l,j

L2(Ω^d). (5.16)

In (5.13) and (5.14), we used the distributive law again and reintroduced the terms we dropped previously. In (5.15), we replaced the ψ_J,iq and ψ_J,jq by the decompositions P_J

lq=1ω_lq,iq and P_J

lq=1ω_lq,jq, respectively. Then, in (5.16), we recombined the product of d one-dimensional L₂-scalar products to one d-dimensional L₂-scalar product. Note that the lower bound with λ_min can be proven in the same way. Now, in combination with (5.9), we know that (5.8) is a norm equivalence with constantsλ^(d),J_max ≤λ^(1),J_max and λ^(d),J_min ≥λ^(1),J_min .

It is almost trivial to prove the sharpness of the estimates, i.e. to show that indeed λ^(d),J_max = λ^(1),J_max and λ^(d),J_min = λ^(1),J_min . Choose the functions u^(1),J_max ∈ V_J^(p), p = 1, . . . , d, for which the

maximum in (5.6) is attained and plug the multivariate function u(x) =

d

Y

p=1

u^(1),J_max(x_p)

into (5.10). This results in an equality instead of an upper bound in (5.13). Theλ^(d),J_min -case can be shown analogously.

We can now easily extend this result to any function u∈H₀¹(Ω^d) by a density argument and λ^(d)_max= lim

J→∞λ^(d),J_max = lim

J→∞λ^(1),J_max =λ⁽¹⁾_max. Again, theλ^(d)_min-case works analogously.

Theorem 5.1 is important to understand why the condition numbers of sparse grid discretiza-tions of the Laplacian decrease with rising dimension. We will discuss this effect in Subsec-tion 5.7.1. The following Theorem also includes a reacSubsec-tion term.

Theorem 5.2. Let λ⁽¹⁾_min and λ⁽¹⁾max be the norm equivalence constants from (5.6) for functions u∈H₀¹(Ω). For u∈H₀¹(Ω^d), αp >0, p= 1, . . . , d,γ ≥0 and

a(u, v) =

d

X

p=1

αp

∂u

∂x_p, ∂v

∂x_p

L2(Ω^d)

+γ(u, v)_L2(Ω^d₎ , (5.17) it holds that

a(u, u)' X

l∈N^d

X^d

p=1

αp2^2lp+ 2γ λ⁽¹⁾_max+λ⁽¹⁾_min

kwlk²_L2(Ω^d₎ for u= X

l∈N^d

wl (5.18)

with w_l ∈W_l,l ∈N^d. The norm equivalence constants λ^(d)_min and λ^(d)max of (5.18)satisfy λ^(d)_min≥ λ⁽¹⁾_min andλ^(d)max≤λ⁽¹⁾max.

Proof. The proof consists of an application of Theorem 5.1 and a simple inequality a(u, u)≤λ⁽¹⁾_maxX

l∈N^d

X^d

p=1

α_p2^2lp

kw_lk²_L2(Ω^d₎+γ(u, v)_L2(Ωd)

=λ⁽¹⁾_maxX

l∈N^d

X^d

p=1

α_p2^2lp

kw_lk²_L2(Ω^d₎+λ⁽¹⁾_max γ λ⁽¹⁾_max

X

l∈N^d

kw_lk²_L2(Ω^d₎

≤λ⁽¹⁾_maxX

l∈N^d

X^d

p=1

α_p2^2lp+ 2γ λ⁽¹⁾_max+λ⁽¹⁾_min

kw_lk²_L2(Ω^d₎.

The lower bound is proven analogously using ^γ

λ⁽¹⁾_min ≥ _λ(1)^2γ min+λ⁽¹⁾max.

5.2 The theory of subspace splittings applied to sparse grids 73

Im Dokument The ANOVA decomposition and generalized sparse grid methods for the high-dimensional backward Kolmogorov equation (Seite 69-81)