Schwarz type algorithms for parabolic problems

Schwarz type algorithms for parabolic problems

Applied Mathematics and Computation 188 (2007) 206–213 www.elsevier.com/locate/amc Schwarz type algorithms for parabolic problems Jianhua Yang School...

169KB Sizes 0 Downloads 45 Views

Applied Mathematics and Computation 188 (2007) 206–213 www.elsevier.com/locate/amc

Schwarz type algorithms for parabolic problems Jianhua Yang School of Information Science and Engineering, Shandong University of Science and Technology, Qingdao 266510, China

Abstract Based on overlapping domain decomposition, we give a parallel subspace correction algorithm for semi-linear parabolic problems. We consider the dependence of convergence rate of this algorithm on parameters of time-step and space-mesh. We give the error estimate, which tell us that the convergence of the approximate solution is independent of the iteration number at each time level. Ó 2006 Elsevier Inc. All rights reserved. Keywords: Parabolic problems; Domain decomposition; Parallel algorithm; Error estimate

1. Introduction Subspace correction method, generally used for construction and analysis of linear iterative methods, is a main technical tool in the analysis of domain decomposition with overlapping. A systematic theory has been developed for elliptic finite element problems in the past few years, see [4,7]. In [2,3] Cai constructed a kind of additive Schwarz algorithms and multiplicative Schwarz method and prove that the convergence rate is smaller than one for parabolic equations. Other domain decomposition methods for parabolic problems can be found in [5,8,9]. In this paper we are interested in solving the parabolic problems using domain decomposition method. We use the Euler backward time discretization and Galerkin approximation in the space variables. At a fixed time level, the resulting equation is equivalent to an elliptic problem which depends on a time step increment and the approximate solution at the last time level. So we can give a parallel algorithm to the parabolic equations at each time level. The key point of this paper is to consider how the convergence and the error estimate depend on diameter of the sub-domain, the special mesh-size, the time step increment and the number of iterations at each time level. The outline of this paper is as follows. In Section 2 we present the time discrete parabolic equations and give a kind of parallel algorithms with Euler backward time discretization. We also give the error estimate result, which tells us that the approximate solution converges after one cycle of iteration at each time level. In Section 3 we give some lemmas. In Section 4 we give the proof of the theorem mentioned in Section 2. E-mail address: [email protected] 0096-3003/$ - see front matter Ó 2006 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2006.09.108

J. Yang / Applied Mathematics and Computation 188 (2007) 206–213

207

Hr(X) denote the general Sobolev spaces with norms kÆk and kÆkr respectively, c and C, with or without subscripts, denote generic positive constants. Their values may be different on different occasions, but are independent of the mesh parameter h, H and time stepping Dt, which will be introduced later. 2. Parallel subspace correction algorithms and convergence results Consider the following parabolic problem: find u(x,t) such that 8  P2 o  ou ou > > < ot  ij¼1 oxj aij ðxÞ oxi ¼ f ðuÞ; ðx; tÞ 2 X  J ; ð1Þ u ¼ 0; ðx; tÞ 2 oX  J ; > > : uðx; 0Þ ¼ u0 ðxÞ; x 2 X; where X is a bounded polygonal domain in R2 with boundary oX and J = (0, T] denote the time interval, aij = aji and there exists a positive constant c such that 2 X 2 X aij ni nj P cjnj2 8n ¼ ðn1 ; n2 ÞT 2 R2 : ð2Þ i¼1

j¼1

The standard variational formulation of the above problem (1) is to find uðtÞ 2 H 10 ðXÞ in t 2 J such that (  ou ; v þ aðu; vÞ ¼ ðf ðuÞ; vÞ; v 2 H 10 ðXÞ; ot v 2 H 10 ðXÞ;

ðuð0Þ; vÞ ¼ ðu0 ; vÞ;

ð3Þ

where aðu; vÞ ¼

Z X 2 X 2 X i¼1

ðf ðuÞ; vÞ ¼

Z

aij

j¼1

ou ov dx; oxi oxj

ð4Þ

f ðuÞv dx:

ð5Þ

X

Let Dt denote the time increment, tn = nDt, un = u(tn), we have  2  un  un1 ou ou þO ¼ Dt : ot ot2 Dt From (3) we can get un 2 H 10 ðXÞ such that  n  u  un1 ; v þ aðu; vÞ ¼ ðf ðun1 Þ; vÞ þ ðqn ; vÞ Dt

8v 2 H 10 ðXÞ;

ð6Þ

ð7Þ

where qn ¼

 2  un  un1 ou ou  þ f ðun Þ  f ðun1 Þ ¼ O Dt : ot ot2 Dt

Then the first equation of (1) can be approximated by  n  u  un1 ; v þ aðu; vÞ ¼ ðf ðun1 Þ; vÞ 8v 2 H 10 ðXÞ: Dt

ð8Þ

Let Aðu; vÞ ¼ ðu; vÞ þ Dtaðu; vÞ8v 2 H 10 ðXÞ; 1=2

kukA ¼ ðAðu; uÞÞ kuka ¼ ðaðu; uÞÞ

1=2

2

2 1=2

¼ ðkuk0 þ Dtkuka Þ

;

:

Then the time discrete problem can be written as AðU n ; V Þ ¼ ðU n1 þ Mtf ðU n1 Þ; V Þ;

8V 2 V:

ð9Þ

208

J. Yang / Applied Mathematics and Computation 188 (2007) 206–213

For a given polygonal region X 2 R2, let Xi be a shape regular finite element triangulation of X and H the N maximal diameter of these X0i s. We sometimes refer Xi as a substructure and fXi gi¼1 as the coarse mesh or H-level triangulation of X. We further divide each substructure Xi into smaller triangles, denoted as sji ; j ¼ 1; . . . We assume that sji form a shape regular finite element triangulation of X and h is the maximal diameter of sji . We callfsji g the fine mesh or h-level triangulation of X. We next define the piecewise linear finite element function spaces over both H-level and h-level triangulations of X VH ¼ fvH j continuous on X; vH jXi linear on Xi ; 8i; vH ¼ 0 on oXg; and Vh ¼ fvh j continuous on X; vh jsj linear on sji ; 8i; vh ¼ 0 on oXg: i

It is obvious that VH  Vh . To obtain a fully discrete problem, we discrete (8) in space by using a Galerkin finite element method and the subspace Vh  H 10 ðXÞ. The Galerkin approximation of (8) reads as follows: Find unh 2 Vh , such that  n  uh  uhn1 ; v þ aðunh ; vÞ ¼ ðf ðuhn1 Þ; vÞ 8v 2 Vh ; ð10Þ Dt ðu0h ; vÞ ¼ ðu0 ; vÞ: ð11Þ 0 0 We extend each subregion Xi to obtain Xi , such that Xi  Xi and there exists a constant a > 0 satisfying: distðoX0i \ X; oXi \ XÞ P aH ; i ¼ 1; 2; . . . ; N : Suppose that oX0i dose not cut through any h-level elements. We make the same constructions for the subregions that meet the boundary except that we cut off the parts that are outside X. To simplify the notations, we N also denote X00 ¼ X. Assume that the decomposition fX0i g1 satisfy the conditions introduced later. Condition (A). For any x 2 X there exists an open domain Di and i 2 {1, 2, . . . , N} such that x 2 Di and Di \ X 2 X0i . Condition (B). For any x 2 X there exists at most N0 subdomain which x can be belong to. It is easy to see that the finite element space Vh can be decomposed into the sum of the coarse mesh function space VH and a number of spaces which are supported only in subregions Xi, i.e. Vh ¼ Vh1 þ Vh2 þ    þ VhN ; where Vhi ¼ Vh \ H 10 ðX0i Þ; i ¼ 1; 2; . . . ; N : Let P hi be the projection from Vh to Vhi with respect to the bilinear form A(Æ, Æ) and P h ¼ P h1 þ P h2 þ    þ P hN . We can obtain the derived equation of (10) as ð12Þ P h uh ¼ g0n1 here g0n1 can be computed without a priori knowledge of the solution uh. It can be seen easily that the operator Ph is usually large and sparse, but not well conditioned. It is wellknown that the efficiency of any iterative methods to be used to solve this linear system of equations depends strongly on the conditioning of the stiffness matrix. The main focus of this paper is to provide a effective preconditioning techniques for solving this equation. From condition (A) we can get an open-overlapping fOi gNi¼1 of X such that Oi \ X 2 P X0i . Under the theory of N unit decomposition (see [1]), there exist smooth function ui (i = 1, 2, . . . , N) satisfying i¼1 ui ¼ 1, P0N 6 ui 6 1, h r 1 h supp(ui)  Oi, and kui kW r;1 6 CH . Let ui ¼ I h ui , here I h : H 0 ðXÞ ! V is a projection, then i¼1 uhi ¼ 1. N Suppose that the overlapping domain decomposition fX0i gi¼1 satisfy condition A and condition B, 0 1 by extending functions in H o ðXi Þ to X by zero, we give the following parallel subspace correction algorithm. Scheme 1. For n P 1 find W n 2 Vh by three steps: 1. Set W n0 ¼ W n1 . 2. Find b e ni;j 2 V hi (j = 1, 2, . . . , m, i = 1, 2, . . . , N) satisfying Aðb e ni;j ; V Þ ¼ ðuhi ðW n1 þ Dtf ðW n1 ÞÞ; V Þ  AðW nj1 ; uhi V Þ 8v 2 Vh :

ð13Þ

J. Yang / Applied Mathematics and Computation 188 (2007) 206–213

209

Let W nj ¼ W nj1 þ

N X

b e nk;j :

k¼1 n

n m,

3. Let W ¼ W x 2 X. Here m denotes the iteration number at each time level. It is clear that the solution of the above scheme is defined uniquely and in Section 4 we can prove the following main result. Theorem 1. Suppose that the solution of (1) is sufficiently smooth and that Wn be the solution of Scheme 1 we have  m ! Dt h 2 kþ1 n n ku  W kA 6 C h þ Dt þ þ ; ð14Þ H2 H where C denotes a generic constant independent of Dt, h and H. 3. Some lemmas In order to prove Theorem 1, in this section we give some auxiliary lemmas. Lemma 1. For all 1 6 i 6 N and 0 6 r 6 1 there exist a constant C0 > 0 such that kðI  I h Þðuhi V ÞkH r 6 C 0

h kV kH r H

8V 2 Vh :

ð15Þ

Proof. Notice the approximation theory of finite element space, we have kðI  I h Þðuhi V ÞkH r 6 Ch2r kD2 ðuhi V ÞkL2 : Since the function uhi and V are all linear, then from the inverse theory of finite element function, we can get kD2 ðuhi V ÞkL2 6 Ckruhi kL1 krV kL2 6 Chr1 H 1 kV kH r ; The proof of Lemma 1 completed.

r ¼ 0; 1:

h

Lemma 2. There exist a constant C1 > 0 such that  ! pffiffiffiffiffi!   N X h Dt   h h þ uk P k W  6 C 1 kV kA kW kA AðV ; W Þ  A V ;   H H k¼1

8V ; W 2 Vh :

Proof. Note that for all V ; W 2 Vh , we have A V;

N X k¼1

! uhk P hk W

¼

N X

Aðuhk V

; P hk W

k¼1

¼ AðV ; W Þ þ

N X k¼1

"   # N 2  X X oV ouhk h ouhk oP hk W Þ þ Dt aij ; P W  aij V ; oxi oxj k oxi oxj k¼1 i;j¼1

AððI  I h Þðuhk V Þ; P hk W Þ

"   # N 2  X X oV ouhk h ouhk oP hk W aij ; P W  aij V ; þ Dt ; oxi oxj k oxi oxj k¼1 i;j¼1

ð16Þ

210

J. Yang / Applied Mathematics and Computation 188 (2007) 206–213

then we have  !  X  N N X   h h uk P k W  6 jAððI  I h Þðuhk V Þ; P hk W Þj AðV ; W Þ  A V ;   k¼1 k¼1    N X 2  h h h X   aij oV ; ouk P h W  aij V ouk ; oP k W  þ Dt k  oxi oxj oxi oxj  k¼1 i;j¼1 6C

N X

kðI  I h Þðuhk V ÞkH 1 kP hk W kA

k¼1

 Dt  krV kL2 kP hk W kL2 þ kV kL2 krP hk W kL2 H pffiffiffiffiffi! h Dt þ kV kA kW kA ; 6C H H þC

where we used Lemma 1. The proof of Lemma 2 completed. h Lemma 3. Suppose that the solution of (1) is sufficiently smooth and that Un be the solution of Galerkin finite element method, then we have kun  unh ks ¼ Oðhrþ1 þ DtÞ;

s ¼ 0; 1:

ð17Þ

This conclusion is well known. 4. Convergence analysis and error estimates In this section, we give the error estimates of the Scheme 1, that is the proof of Theorem 1. Proof of Theorem 1. In order to get the error estimates, first we give the error equation. Actually the parallel subspace correction algorithm is a iterative method to solve the following equation with a initial approximation Wn1 b n ; V Þ ¼ ðDtf ðW n1 Þ þ W n1 ; V Þ Að W

8V 2 Vh ;

ð18Þ

b n  U n and en = Wn  Un, we have Let E ¼ W n

AðEn ; V Þ ¼ ðen1 ; V Þ þ Dtðf ðW n1 Þ  f ðU n1 Þ; V Þ 8V 2 V h

ð19Þ

or equivalently, b n ; V Þ þ ðen1 ; V Þ þ Dtðf ðW n1 Þ  f ðU n1 ÞÞ; V Þ Aðen ; V Þ ¼ AðW n  W

8V 2 V h :

b n. From the error equation, we can see that the key is to estimate W n  W P n n From the iterative equation of the parallel algorithm: W j ¼ W j1 þ Ni¼0 b e ni;j , we have that ! N X b n ; V Þ ¼ AðW n  W b n; V Þ þ A b AðW n  W 8V 2 Vh : en ; V j

j1

i;j

ð20Þ

ð21Þ

i¼0

On the other hand we have A

N X i¼0

! be ni;j ; V

¼

N h i X uhi ðW n1 þ Dtf ðW n1 ÞÞ; P hi V Þ  AðW nj1 ; uhi P hi V Þ i¼0

N h i X b n ; uh P h V Þ  AðW n ; uh P h V Þ : ¼ Að W i i j1 i i i¼0

ð22Þ

J. Yang / Applied Mathematics and Computation 188 (2007) 206–213

Let enj ¼ W nj  U n , then bnWn ; Aðenj ; V Þ ¼ Aðenj1 ; V Þ þ A W j1

N X

211

! 8V 2 Vh :

uhi P hi V

ð23Þ

i¼0

So we have Aðenj

n

 E ;VÞ ¼

Aðenj1

n

n

 E ;VÞ þ A E 

enj1 ;

N X

! uhi P hi V

¼A

enj1

n

E ; I 

i¼0

From all of the above, we can get Aðenj

E

n

; enj

n

E Þ¼A

enj1

n

E ; I 

N X

! ! uhi P hi

V

:

ð24Þ

i¼0

!

! uhi P hi

N X

ðenj

n

E Þ :

ð25Þ

i¼0

Combining these inequalities with the inequality of Lemma 2, we have  2  h Dt 2 2 kenj  En kA 6 C 1 þ kenj1  En kA ; 1 6 j 6 m; H2 H2

ð26Þ

we derive that n

ke 

2 E n kA





h2 Dt þ 2 2 H H

6 C1

 m

2

ken0  En kA :

ð27Þ

Noticing 2

2

2

ken0  En kA ¼ kW n1  un1 þ uhn1  unh  En kA ¼ ken1  En þ U n1  U n kA ; h so we have that n



n



ke  E kA 6 C 1

h2 Dt þ H2 H2

 m2

ke

n1

oU k dt :  E kA þ k ot A tn1 n

Z

tn

ð28Þ

Hence

 2  m

Z tn oU h Dt 2 n1 n dt : þ ke  E k þ ken kA 6 kEn kA þ C 1 A ot 2 2 H H tn1 A

ð29Þ

Now we estimate the terms on the right-hand side of (29). Letting v = En in (19) we have 2

kEn kA ¼ ðen1 ; En Þ þ Dtðf ðW n1 Þ  f ðU n1 Þ; En Þ ¼ ðen1 ; En Þ þ Dtðf 0  en1 ; En Þ 6 ð1 þ CDtÞken1 kkEn k; 0

where f denotes the first order derivative of f at a point between W so we have

ð30Þ n1

n1

and U

0

,f Æe

kEn kA 6 ð1 þ CDtÞken1 k:

n1

= f(W

n1

)  f(U

n1

).

ð31Þ

From (19), we have that ðEn  en1 ; vÞ þ DtaðEn  en1 ; vÞ ¼ Dtðf 0  en1 ; vÞ  Dtaðen1 ; vÞ: n

Let v = E  e

n1

ð32Þ

, similarly to (31) we can prove that

kEn  en1 kA 6 Dt1=2 ken1 ka þ CDtken1 k: From (29), (31) and (33) we have that  2 m=2  2 m=2 h Dt h Dt 1=2 n1 þ Dt ke k þ CDt þ : ken kA 6 ð1 þ CDtÞken1 k þ C a H2 H2 H2 H2

ð33Þ

ð34Þ

212

J. Yang / Applied Mathematics and Computation 188 (2007) 206–213

By 2ab 6 a2 + b2 we have that !2  2 m=2 h Dt 1=2 n1 n1 þ ðDtÞ ke ka ke k þ C H2 H2  2 m=2  2 m h Dt Dt 1=2 2 2 h n1 2 n1 n1 ¼ ke k þ 2 þ ðDtÞ ke ka ke k þ C þ Dtken1 ka H2 H2 H2 H2  2 m  2 m h Dt Dt 2 h n1 6 ken1 k2 þ Dtken1 k2a þ C 2 þ ke k þ C þ Dtken1 k2a H2 H2 H2 H2  2 m ! Dt 2 2 h 6 1þC þ ken1 kA : H2 H2 It is clear that there exists a constant C5 such that  2 m !  2 m Dt h Dt 2 h þ þ 6 1þC 6 1 þ 2C 5 H2 H2 H2 H2

 1 þ C5

h2 Dt þ 2 2 H H

ð35Þ

m !2

therefore 1þC

2



h2 Dt þ 2 2 H H

m !1=2

 6 1 þ C5

h2 Dt þ 2 2 H H

m

so we have ke

n1



h2 Dt kþC þ H2 H2

m=2 ðDtÞ

1=2

ke

n1

 ka 6

1 þ C5

h2 Dt þ H2 H2

m ! ken1 kA :

Using (34) and (35) we have that  2 m !  2 m=2 h Dt h Dt n n1 ke kA 6 1 þ C 5 ke kA þ CDt þ þ ; H2 H2 H2 H2

ð36Þ

where the constant C is independent of Dt, h and H. As discussed in [6], we have that  m ! Dt h2 2 kþ1 n n ku  W kA 6 C h þ Dt þ þ : H2 H2

ð37Þ

We complete the proof of Theorem 1. h Remark 1. Lemma 3 and (37) tell us that the convergence rate of Scheme 1 at each time-level is C



Dt H2

12 2 þ Hh 2 .

Acknowledgement The author gives her thanks to Professor Yang Danping for his suggestion and help. References [1] [2] [3] [4] [5]

R.A. Adamas, Sobolev Space, Academic Press, New York, 1975. X.C. Cai, Multiplicative Schwarz methods for parabolic problem, SIAM J. Sci. Comput. 15 (1994) 587–603. X.C. Cai, Additive Schwarz algorithms for parabolic convection diffusion equations, Numer. Math. 60 (1991) 41–61. J. Xu, Iterative methods by space decomposition and subspace correction, SIAM Rev. 34 (4) (1992) 581–613. C. Dawson, Q. Du, T.F. Dupont, A finite difference domain decomposition algorithm for numerical solution of the heat equation. Tech. Rep., 89-09, Univ. of Chicago, 1989.

J. Yang / Applied Mathematics and Computation 188 (2007) 206–213

213

[6] J. Yang, D. Yang, Additive Schwarz methods for parabolic problems, Appl. Math. Comput. 163 (2005) 17–28. [7] J. Xu, Iterative methods by SPD and small subspace solvers for nonsymmetric or indefinite problems, in: D.E. Keyes, T.F. Chan, G.A. Meurant, J.S. Scroggs, R.G. Voigt (Eds.), fifth International Symposium on Domain Decomposition Methods for Partial Differential Equations, SIAM, Philadelphia, 1992, pp. 106–118. [8] R.E. Ewing, R.D. Lazarov, J.E. Pasciak, P.S. Vassilevski, Finite element methods for parabolic problems with time steps variable in space. Tech. Rep., 1989-05, Inst. for Sci. Comp., Univ. of Wyoming, 1989. [9] Yu.A. Kuznetsov, Domain decomposition methods for time dependent problems, Preprint, 1989.