Fashion & Beauty

15 pages
5 views

A local discontinuous Galerkin method for the second-order wave equation

of 15
All materials on our website are shared by users. If you have any questions about copyright issues, please report us to resolve them. We are always happy to assist you.
Share
Description
A local discontinuous Galerkin method for the second-order wave equation
Transcript
  A local discontinuous Galerkin method for the second-order wave equation Mahboub Baccouch ⇑ Department of Mathematics, University of Nebraska at Omaha, Omaha, NE 68182, United States a r t i c l e i n f o  Article history: Received 28 February 2011Received in revised form 24 August 2011Accepted 21 October 2011Available online 2 November 2011 Keywords: Local discontinuous Galerkin methodSecond-order wave equationSuperconvergence  A posteriori  error estimation a b s t r a c t In this paper we present new superconvergence results for the local discontinuous Galerkin (LDG)method applied to the second-order scalar wave equation in one space dimension. Numerical experi-ments show  O ð h  p þ 1 ÞL 2 convergence rate for the LDG solution and  O ( h  p +2 ) superconvergent solutions atRadau points. More precisely, a local error analysis reveals that, at a fixed time  t  , the leading terms of the discretization errors for the solution and its derivative using  p -degree polynomial approximationsare proportional to the (  p  + 1)-degree right Radau and (  p  + 1)-degree left Radau polynomials, respectively.Thus, the  p -degree LDG solution is  O ( h  p +2 ) superconvergent at the roots of the (  p  + 1)-degree right Radaupolynomial and the derivative of the LDG solution is  O ( h  p +2 ) superconvergent at the roots of the (  p  + 1)-degree left Radau polynomial. These results are used to construct simple, efficient, and asymptoticallycorrect  a posteriori  error estimates in regions where solutions are smooth. Finally, we present severalnumerical examples to validate the superconvergence results and the asymptotic exactness of our  a pos-teriori  errors estimates under mesh refinement.   2011 Elsevier B.V. All rights reserved. 1. Introduction The wave equation is important in the description of many trav-eling wave phenomena such as the vibrations of an elastic bar, thesound and light waves, and so on. The numerical solution of thewave equation is of fundamental importance to the simulation of time-dependent sound, light, radiation of electromagnetic, or elas-tic waves. For such wave phenomena the second-order wave equa-tion often serves as a model problem. Standard numerical methodssuch as finite difference, continuous Galerkin finite element, and fi-nite volume methods have poor stability properties when appliedto hyperbolic problems. In contrast, the discontinuous Galerkin(DG) method is known to have good stability properties when ap-plied to hyperbolic problems.The DG finite element method was first used to solve the neu-tron equation [23]. Cockburn and Shu [15] extended the method to solve first-order hyperbolic partial differential equations of con-servation laws. They also developed the local discontinuous Galer-kin (LDG) method for convection–diffusion problems [16]. Consult[14] and the referencescitedthereinfor a detaileddiscussionof thehistory of DG method and a list of important citations on the DGmethod and its applications.DG methods allow discontinuous bases, which simplify both h -refinement (mesh refinement and coarsening) and  p -refine-ment (method order variation). However, for DG methods tobe used in an adaptive framework one needs  a posteriori  errorestimates to guide adaptivity and stop the refinement process.The solution space consists of piecewise continuous polynomialfunctions relative to a structured or unstructured mesh. As such,it can sharply capture solution discontinuities relative to thecomputational mesh. It maintains local conservation on an ele-mental basis. The success of the DG method is due to the follow-ing properties: ( i ) does not require continuity across elementboundaries, ( ii ) is locally conservative, ( iii ) is well suited to solveproblems on locally refined meshes with hanging nodes, ( i v  )exhibits strong superconvergence that can be used to estimatethe discretization error, ( v  ) has a simple communication patternbetween elements with a common face that makes it useful forparallel computation and ( vi ) it can handle problems with com-plex geometries to high order. Adjerid and Baccouch [1,2,7]investigated the superconvergence properties of the DG methodapplied to steady and transient hyperbolic partial differentialequations on structured and unstructured triangular meshes.They used the superconvergence results to construct simple, effi-cient, and asymptotically correct  a posteriori  error estimates bysolving local hyperbolic problems with no boundary conditionson triangular meshes.The LDG finite element method is an extension of the DG meth-od aimed at solving partial differential equations containing higherthan first order spatial derivatives. The LDG method possesses sev-eral properties which make it popular for practical computations.As it is well known, the success of the LDG method is due to thefollowing properties: ( i ) is locally conservative, ( ii ) allows a veryefficient parallelization and ( iii ) is suitable for  h -,  p -, and  hp -refine- 0045-7825/$ - see front matter    2011 Elsevier B.V. All rights reserved.doi:10.1016/j.cma.2011.10.012 ⇑ Tel.: +1 402 554 4016; fax: +1 402 554 2975. E-mail address:  mbaccouch@unomaha.eduComput. Methods Appl. Mech. Engrg. 209–212 (2012) 129–143 Contents lists available at SciVerse ScienceDirect Comput. Methods Appl. Mech. Engrg. journal homepage: www.elsevier.com/locate/cma  ments. Consult [10] and the references citedtherein for more infor-mation on LDG methods.The main motivation for the LDG method proposed in this papersrcinates from the LDG techniques which have been developed forconvection–diffusion equations. The LDG method for solving con-vection–diffusion problems was first introduced by Cockburn andShu in [16]. They further studied the stability and error estimatesfor the LDG method. Castillo et al. [9] presented the first  a priori  er-ror analysis for the LDG method for a model elliptic problem. Cock-burn et al. [13] presented a superconvergence result for the LDGmethod for a model elliptic problem on Cartesian grids. Pointwisesuperconvergence properties of LDG methods for one-dimensionaldiffusion problems were investigated by Castillo [8], Adjerid andIssaev [5], Celiker and Cockburn [11], Zhang et al. [26], and Cheng and Shu [12].  A posteriori  error estimates lie in the heart of every adaptive fi-nite element algorithm for differential equations. They are used toassess the quality of numerical solutions and guide the adaptiveenrichment process where elements having high errors are en-riched by  h -refinement and/or  p -refinement while elements withsmall errors are  h -and/or  p -coarsened. Furthermore, error esti-mates are used to stop the adaptive refinement process. An  a pos-teriori  error analysis for convection and convection-dominatedconvection–diffusion problems was presented in [22,19]. Several a posteriori  DG error estimates are known for diffusive [20,24]problems. Adjerid and Baccouch [3] investigated the global conver-gence of the implicit residual-based  a posteriori  error estimates of Adjerid et al. [4] for a scalar linear hyperbolic problem in onedimension. They proved that, for smooth solutions, these  a posteri-ori  error estimates at a fixed time  t   converge to the true spatial er-ror in the  L 2 -norm under mesh refinement.The LDG method for transient convection–diffusion problemswere investigated by Adjerid and Klauser [6]. The authors pre-sented a numerical study of superconvergence properties for theLDG method of Cockburn and Shu [16] applied to time-dependentconvection–diffusion problems in one dimension. In particular,they observed that  p -degree piecewise polynomial discontinuousfinite element solutions of convection-dominated problems are O ( h  p +2 ) superconvergent at the shifted roots of the (  p  + 1)-degreeRadau polynomial on each element. For diffusion-dominated prob-lems, they observed that the LDG errors in  u  and  u  x  do not exhibitany superconvergence, while they observed that the derivative of the LDG solution is  O ( h  p +2 ) superconvergent at the roots of thederivative of Radau polynomial of degree  p  + 1. Furthermore, theyassumed that the leading term of the discretization error for thederivative  u  x  using  p -degree polynomial approximations is propor-tional to the (  p  + 1)-degree right Radau polynomials. They usedthese results to construct asymptotically exact  a posteriori  finiteelement error estimates.The second-order wave equation can be reformulated as afirst-order hyperbolic system of equations. Many DG schemes forsolving first-order hyperbolic systems are available [20,24]. For in-stance, Cockburn and Shu [15] used a DGFEM in space combinedwith a Runge–Kutta scheme in time to discretize hyperbolic con-servation laws. Hesthaven and Warburton [18] used the same ap-proach to implement high-order methods for Maxwell’s equations. Johnson [21] implemented the DG method for second-order hyper-bolic problems in first-order hyperbolic form. Consult [21,17,25]and the references cited therein for more details and for othernumerical schemes developed to solve the wave equation usingDG methods.In this paper we analyze the spatial discretization errors associ-ated with solutions of one-dimensional second-order wave equa-tion by the LDG method in space. To the best knowledge of theauthor, this is the first superconvergent LDG method for the sec-ond-order wave equation. Numerical experiments show O ð h  p þ 1 ÞL 2 convergence rate for the LDG solution and  O ( h  p +2 ) super-convergent solutions at Radau points. More precisely, a local erroranalysis reveals that the leading terms of the discretization errorsfor the solution  u  and its derivative  u  x  using  p -degree polynomialapproximations are proportional to the (  p  + 1)-degree right Radauand (  p  + 1)-degree left Radau polynomials, respectively. Thus, the  p -degree discontinuous finite element solutionis  O ( h  p +2 ) supercon-vergent at the roots of the (  p  + 1)-degree right Radau polynomialand the derivative of the LDG solution is  O ( h  p +2 ) superconvergentat the roots of the (  p  + 1)-degree left Radau polynomial. The globalsuperconvergence results are numerically confirmed. We applythese superconvergence results to construct a superconvergence-based  a posteriori  LDG error estimates for the second-order waveequation. These  a posteriori  error estimates are simple, efficientand asymptotically correct. The leading terms of the discretizationerrors for  u  and  u  x  are estimated by solving local problems with noboundary conditions on each element of the mesh.The paper is organized as follows. In Section 2, we present themodel problem and define the LDG method for the wave equation.Details related to the implementation of the method are also de-scribed in Section 2. In Section 3, we present the local error analy- sis and new superconvergence results. In Section 4, we show our  a posteriori  error estimation procedure for the second-order waveequation. In Section 5, we present numerical results and show thatthe global error estimate converges to the true error under meshrefinement in  L 2 . Concluding remarks are given in Section 6. 2. The LDG method for the wave equation Let  u (  x , t  ) denote a smooth function on [  1,1] and consider thefollowing linear second-order scalar wave equation in one spacevariable u tt     c  2 u  xx  ¼  f  ð  x ; t  Þ ;  x  2 ð 1 ; 1 Þ ;  t   2 ð 0 ; T   ;  ð 2 : 1a Þ subject to the initial and Dirichlet boundary conditions u ð  x ; 0 Þ ¼  g  ð  x Þ ;  u t  ð  x ; 0 Þ ¼  h ð  x Þ ;  x  2 ½ 1 ; 1  ;  ð 2 : 1b Þ u ð 1 ; t  Þ ¼  u l ð t  Þ ;  u ð 1 ; t  Þ ¼  u r  ð t  Þ ;  t   2 ½ 0 ; T   :  ð 2 : 1c Þ We will also investigate (2.1a) subject to the initial conditions(2.1b) and the mixed boundary conditions u ð 1 ; t  Þ ¼  u l ð t  Þ ;  u  x ð 1 ; t  Þ ¼  u r  ð t  Þ ;  t   2 ½ 0 ; T   :  ð 2 : 1d Þ Typically  u (  x , t  ) represents the displacement at point  x  and at time  t  of an elastic bar from its equilibrium position in the  y -direction. Theconstant  c   is the wave speed and [0, T  ] is a finite time interval. In ouranalysis we assume that the constant  c   > 0 and we select the initialand boundary conditions and the (known) source,  f  (  x , t  ), such thatthe exact solution,  u (  x , t  ), is a smooth function. We also assume thatthe initial conditions are consistent with the boundary data. In order to construct the LDG formulation, we introduce an aux-iliary variable  q (  x , t  ) =  u  x (  x , t  ) and rewrite our model problem (2.1)as a first-order system u tt     c  2 q  x  ¼  f  ð  x ; t  Þ ;  x  2 ð 1 ; 1 Þ ;  t   2 ð 0 ; T   ;  ð 2 : 3a Þ q   u  x  ¼  0 ;  x  2 ð 1 ; 1 Þ ;  t   2 ð 0 ; T   :  ð 2 : 3b Þ  2.1. Discretization in space InordertoobtaintheweakLDGformulation,wepartition[  1,1]into a quasi-uniform mesh, {  1 =  x 0  <  x 1  <  x 2  <  <  x N   1  <  x N   = 1},having  N   subintervals I  i  = [  x i  1 ,  x i ],  i  = 1, . . . , N  . Throughoutthis man-uscript  h i  =  x i   x i  1  and  h  = max 1 6 i 6 N   h i  such that  hh i <  K  , where K   > 1 independent of   h . 130  M. Baccouch/Comput. Methods Appl. Mech. Engrg. 209–212 (2012) 129–143  Let us multiply (2.3a) and (2.3b) by arbitrary smooth test func-tions  v   and  w , respectively, integrate over an arbitrary element  I  i ,and use integration by parts to write Z  I  i u tt  v   dx þ  c  2 Z  I  i q v  0 dx    c  2 q ð  x i ; t  Þ v  ð  x i Þ þ  c  2 q ð  x i  1 ; t  Þ v  ð  x i  1 Þ¼ Z  I  i  f  v   dx ;  8  v   2 L 2 ;  ð 2 : 4a Þ Z  I  i qwdx  þ Z  I  i uw 0 dx    u ð  x i ; t  Þ w ð  x i Þþ  u ð  x i  1 Þ w ð  x i  1 ; t  Þ¼  0 ;  8  w  2 L 2 :  ð 2 : 4b Þ Next, we construct the piecewise-polynomial space  V   ph  as the spaceof polynomials of the degree up to  p  in each cell  I  i , i.e. V   ph  ¼  v   :  v  j I  i 2 P   p ;  for  x  2  I  i ;  i  ¼  1 ; . . . ; N  no ;  ð 2 : 5a Þ where  P   p  is the space of polynomials of degree at most  p . Next, we approximate the exact solutions  u  and  q  by piecewisepolynomials  u h  2 V   ph  and  q h  2 V   ph , respectively, whose restriction to I  i  are in P   p . Note that the functions  u h  and  q h  are piecewise polyno-mials not necessarily continuous at the nodes of the subintervals  I  i .The semi-discrete LDG methodconsistsof finding u h  and  q h  2 V   ph such that, for all  v   and  w  2 V   ph , Z  I  i ð u h Þ tt  v   dx þ  c  2 Z  I  i q h v  0 dx    c  2 ^ q h ð  x i ; t  Þ v  ð  x  i  Þþ  c  2 ^ q h ð  x i  1 ; t  Þ v  ð  x þ i  1 Þ¼ Z  I  i  f  v   dx ;  ð 2 : 6a Þ Z  I  i q h wdx þ Z  I  i u h w 0 dx  ^ u h ð  x i ; t  Þ w ð  x  i  Þþ ^ u h ð  x i  1 ; t  Þ w ð  x þ i  1 Þ¼ 0  ð 2 : 6b Þ for all  i  = 1, . . . , N  , subject to the initial conditions  u h ð  x ; 0 Þ ¼  P  h u 0 ð  x Þ2 V   ph  and  @  u h @  t   ð  x ; 0 Þ ¼  P  h u 1 ð  x Þ 2 V   ph  obtained by using the local  L 2 pro- jection on  V   ph ,  i.e. , Z  I  i P  h u 0 v  dx  ¼ Z  I  i  g  v  dx ;  8  v   2 V   ph ; Z  I  i P  h u 1 v  dx  ¼ Z  I  i h v  dx ;  8  v   2 V   ph :  ð 2 : 6c Þ We will also compute  u h ð  x ; 0 Þ 2 V   ph  and  @  u h @  t   ð  x ; 0 Þ 2 V   ph  by interpolat-ing  g  (  x ) and  h (  x ) as u h ð  x ; 0 Þ ¼  p  g  ð  x Þ ; @  u h @  t   ð  x ; 0 Þ ¼  p h ð  x Þ ;  x  2  I  i ;  ð 2 : 6d Þ where  p w  is the  p -degree polynomial that interpolates  w  at theroots of (  p  + 1)-degree right Radau polynomial. The numerical fluxes  ^ u h  and  ^ q h  are the discrete approximationsto the traces of   u  and  q , respectively, at all nodes. They should bedesigned based on different guiding principles for the partial dif-ferential equation to ensure stability.We use the usual notations to denote the average, { v  }, and jump, [ v  ], values of a function  v   at  x i f v  gð  x i Þ ¼  12 ð v  þ ð  x i Þþ  v   ð  x i ÞÞ ;  ½ v  ð  x i Þ ¼  v  þ ð  x i Þ  v   ð  x i Þ ;  ð 2 : 6e Þ where  v  + and  v   are the values of   v   at  x i  from the right cell  I  i +1 , andfrom the left cell  I  i , respectively, v   ð  x i Þ ¼  lim s ! 0   v  ð  x i  þ  s Þ ;  v  þ ð  x i Þ ¼  lim s ! 0 þ v  ð  x i  þ  s Þ :  ð 2 : 6f  Þ In order to complete the definition of the LDG method we need toselect the numerical fluxes  ^ u h  and  ^ q h  on the boundaries of   I  i . In thispaper, we take the following numerical fluxes 1. The (conservative and consistent) numerical fluxes  ^ u h  and  ^ q h  forthe Dirichlet boundary conditions (2.1c) are defined by ^ u h ð  x i ; t  Þ ¼ u l ð t  Þ ;  i  ¼  0 ; u  h ð  x i ; t  Þ ;  i  ¼  1 ; . . . ; N    1 ; u r  ð t  Þ ;  i  ¼  N  ; 8><>: ; ^ q h ð  x i ; t  Þ ¼ q þ h ð  x i ; t  Þ ;  i  ¼  0 ; . . . ; N    1 ; q  h ð  x N  ; t  Þ þ  d ð u  h ð  x N  ; t  Þ   u r  ð t  ÞÞ ;   ð 2 : 6g Þ where the stabilization parameter,  d , for our LDG method is givenby  d  ¼  ph i . 2. The numerical fluxes for the mixed boundary conditions (2.1d)are chosen as ^ u h ð  x i ; t  Þ ¼ u l ð t  Þ ;  i  ¼  0 ; u  h ð  x i ; t  Þ ;  i  ¼  1 ; . . . ; N  :   ; ^ q h ð  x i ; t  Þ ¼  q þ h ð  x i ; t  Þ ;  i  ¼  0 ; . . . ; N    1 ; u r  ð t  Þ ;  i  ¼  N  :   ð 2 : 6h Þ In our analysis we need Jacobi polynomials defined by Rodriguesformula P  a ; b n  ð n Þ ¼ ð 1 Þ n 2 n n ! ð 1   n Þ  a ð 1 þ  n Þ  b  d n d n n  ð 1   n Þ a þ n ð 1 þ  n Þ b þ n hi ; a ; b  >   1 ;   1 6 n 6 1 ;  ð 2 : 7a Þ where  L n ð n Þ ¼  P  0 ; 0 n  ð n Þ ;  1 6 n 6 1, is the  n th-degree Legendre poly-nomial and satisfies the orthogonality relation Z   1  1 L n ð n Þ L m ð n Þ d n  ¼  2 d nm 2 n  þ 1 ;  n ; m P 0 ;  ð 2 : 7b Þ where  d nm  is the Kronecker symbol equal to 1 if   n  =  m  and 0 other-wise. The (  p  + 1)th-degree right Radau polynomial R þ  p þ 1 ð n Þ ¼  L  p þ 1 ð n Þ  L  p ð n Þ ;   1 6 n 6 1  ð 2 : 7c Þ has  p  + 1 real distinct roots,   1  <  n þ 0  <  n þ 1  <    <  n þ  p  1  <  n þ  p  ¼  1. The(  p  + 1)th-degree left Radau polynomial R   p þ 1 ð n Þ ¼  L  p þ 1 ð n Þþ  L  p ð n Þ ;   1 6 n 6 1  ð 2 : 7d Þ has  p  + 1 real distinct roots,   1  ¼  n  0  <  n  1  <  n  2  <    <  n   p  1 <  n   p  <  1. In this paper, we define the  L 2 inner product of two integrablefunctions  u  and  v   on  I  i  and [  1,1] as ð u ; v  Þ i  ¼ Z   x i  x i  1 u ð  x Þ v  ð  x Þ dx ;  ð u ; v  Þ ¼ Z   1  1 u ð n Þ v  ð n Þ d n  ð 2 : 8a Þ and the subsequent induced norms are k u k 2 i  ¼ ð u ; u Þ i ;  k u k 2 ¼ ð u ; u Þ :  ð 2 : 8b Þ  2.2. Discretization in time Expressing the LDG solutions  u h  and  q h  as a linear combinationof orthogonal basis u h ð  x ; t  Þ ¼ X  pk ¼ 0 c  k ð t  Þ L k ð  x Þ ;  q h ð  x ; t  Þ ¼ X 2  p þ 1 k ¼  p þ 1 c  k ð t  Þ L k   p  1 ð  x Þ ;  x  2  I  i ð 2 : 9 Þ and testing against functions  v   =  w  =  L  j (  x ) in (2.6a) and (2.6b), thediscretization of  (2.3) in space by the LDG method (2.6) leads to the following equations M d 2 C  ! u dt  2  ¼  F  !ð t  ;  C  ! u ;  C  ! q Þ ;  t   2 ð 0 ; T   ;  ð 2 : 10a Þ G !ð t  ;  C  ! u ;  C  ! q Þ ¼ ~ 0 ;  t   2 ½ 0 ; T   ;  ð 2 : 10b Þ M. Baccouch/Comput. Methods Appl. Mech. Engrg. 209–212 (2012) 129–143  131  subject to the initial conditions C  ! u ð 0 Þ ¼  ~ a u ;  dC  ! u dt   ð 0 Þ ¼ ~ b u ;  ð 2 : 10c Þ where  C  ! u ð t  Þ ¼ ½ c  0 ð t  Þ ; . . . ; c   p ð t  Þ t  ;  C  ! q ð t  Þ ¼ ½ c   p þ 1 ð t  Þ ; . . . ; c  2  p þ 1 ð t  Þ t  , and M   = ( m kj ) is the mass matrix with elements  m kj  ¼ R  I  i L k L  j dx . The vec-tors  F  !ð t  ;  C  ! u ;  C  ! q Þ  and  G !ð t  ;  C  ! u ;  C  ! q Þ  are given by F  !ð t  ;  C  ! u ;  C  ! q Þ ¼ Z  I  i  f  ~ / dx    c  2 Z  I  i q h ~ / 0 dx þ  c  2 ^ q h ð  x i ; t  Þ ~ /  ð  x i Þ   ^ q h ð  x i  1 ; t  Þ ~ / þ ð  x i  1 Þ  ;  ð 2 : 10d Þ G !ð t  ;  C  ! u ;  C  ! q Þ ¼ Z  I  i q h ~ w dx þ Z  I  i u h ~ w 0 dx    ^ u h ð  x i ; t  Þ ~ w  ð  x i Þþ  ^ u h ð  x i  1 ; t  Þ ~ w þ ð  x i  1 Þ ;  ð 2 : 10e Þ where  ~ / ð  x Þ ¼ ½ L 0 ð  x Þ ; . . . ; L  p ð  x Þ t  and  ~ w ð  x Þ ¼ ½ L 0 ð  x Þ ; . . . ; L  p ð  x Þ t  . The initial values  ~ a u  and ~ b u  are obtained from the initial condi-tions  u (  x ,0) and  u t  (  x ,0), respectively, using either the  L 2 projectionor by interpolating  g  (  x ) and  h (  x ) at the roots of (  p  + 1)-degree rightRadau polynomial as described in (2.6c) and (2.6d).To discretize (2.10a) in time, we transform (2.10a) into a system of first-order ordinary differential equations. To do this, we definenew variables  Y  ! 1  and  Y  ! 2  as follows Y  ! 1  ¼  C  ! u ;  Y  ! 2  ¼  dC  ! u dt   ;  Y  ! 3  ¼  C  ! q and write (2.10) as dY  ! 1 dt   ¼  Y  ! 2 ;  t   2 ð 0 ; T   ;  ð 2 : 11a Þ dY  ! 2 dt   ¼  M   1 F  !ð t  ;  Y  ! 1 ;  Y  ! 3 Þ ;  t   2 ð 0 ; T   ;  ð 2 : 11b Þ G !ð t  ;  Y  ! 1 ;  Y  ! 3 Þ ¼ ~ 0 ;  t   2 ½ 0 ; T   ;  ð 2 : 11c Þ subject to the initial conditions Y  ! 1 ð 0 Þ ¼  ~ a u ;  Y  ! 2 ð 0 Þ ¼ ~ b u :  ð 2 : 11d Þ Using the notation S  !ð t  ;  Y  ! ;  Y  ! 3 Þ ¼  Y  ! 2 M   1 F  !ð t  ;  Y  ! 1 ;  Y  ! 3 Þ "# ;  Y  !ð t  Þ ¼  Y  ! 1 ð t  Þ Y  ! 2 ð t  Þ "# ; the system (2.11) can be solved for  c  k ,  k  = 0, . . . ,2  p  + 1 using  e.g. , theclassical fourth order Runge–Kutta method. Let D t   denote the timestep and set  t  n  =  n D t  . Then the Runge–Kutta method consists infinding approximations Y  ! ð n Þ ¼  Y  ! ð n Þ 1 Y  ! ð n Þ 2 "#  to  Y  !ð t  n Þ ¼  Y  ! 1 ð t  n Þ Y  ! 2 ð t  n Þ "# ; using the following algorithm:  Y  ! ð 0 Þ ¼  Y  !ð 0 Þ ¼ ½ ~ a u ;~ b u  t  is known. For n  = 0,1,2, . . . , Compute  Y  ! ð n Þ 3  from  G !ð t  n ;  Y  ! ð n Þ ;  Y  ! ð n Þ 3  Þ ¼ ~ 0 and set  K  ! 1  ¼  S  !ð t  n ;  Y  ! ð n Þ ;  Y  ! ð n Þ 3  Þ ; Compute  Y  ! ð n Þ 3  from  G !ð t  n  þ D t  2  ;  Y  ! ð n Þ þ D t  2  K  ! 1 ;  Y  ! ð n Þ 3  Þ ¼ ~ 0andset  K  ! 2  ¼  S  ! t  n  þ D t  2  ;  Y  ! ð n Þ þ D t  2  K  ! 1 ;  Y  ! ð n Þ 3  ; Compute  Y  ! ð n Þ 3  from  G !ð t  n  þ D t  2  ;  Y  ! ð n Þ þ D t  2  K  ! 2 ;  Y  ! ð n Þ 3  Þ ¼ ~ 0and set  K  ! 3  ¼  S  ! t  n  þ D t  2  ;  Y  ! ð n Þ þ D t  2  K  ! 2 ;  Y  ! ð n Þ 3  ; Compute  Y  ! ð n Þ 3  from  G !ð t  n  þ D t  ;  Y  ! ð n Þ þ D tK  ! 3 ;  Y  ! ð n Þ 3  Þ ¼ ~ 0and set  K  ! 4  ¼  S  ! t  n  þ D t  ;  Y  ! ð n Þ þ D tK  ! 3 ;  Y  ! ð n Þ 3  ; Y  ! ð n þ 1 Þ ¼  Y  ! ð n Þ þ D t  6  K  ! 1  þ 2 K  ! 2  þ 2 K  ! 3  þ  K  ! 4  ;  n  ¼ 0 ; 1 ; 2 ; . . . ð 2 : 12 Þ Let the local LDG errors  e u  and  e q  be defined as e u ð  x ; t  Þ ¼  u ð  x ; t  Þ  u h ð  x ; t  Þ ;  e q ð  x ; t  Þ ¼  q ð  x ; t  Þ  q h ð  x ; t  Þ ; where  u h  and  q h  are the LDG solutions on  I  i  = [  x i  1 ,  x i ]. We subtract (2.6a) from (2.4a) with  v   2 V   ph  and (2.6b) from(2.4b) with  w  2 V   ph  to obtain the LDG orthogonality conditions forthe local errors  e u  and  e q  on  I  i Z  I  i ð e u Þ tt  v   dx  þ  c  2 Z  I  i e q v  0 dx    c  2 ^ e q ð  x i ; t  Þ v   ð  x i Þþ c  2 ^ e q ð  x i  1 ; t  Þ v  þ ð  x i  1 Þ ¼  0 ;  ð 2 : 13a Þ Z  I  i e q wdx  þ Z  I  i e u w 0 dx   ^ e u ð  x i ; t  Þ w  ð  x i Þþ ^ e u ð  x i  1 ; t  Þ w þ ð  x i  1 Þ ¼  0 :  ð 2 : 13b Þ Using the mapping of   I  i  = [  x i  1 ,  x i ] onto the canonical element [  1,1]defined by  x ð n ; h Þ ¼  x i  þ  x i  1 2  þ  h 2 n ;  ð 2 : 14 Þ and letting ~ e u ð n ; h ; t  Þ ¼  u ð  x ð n ; h Þ ; t  Þ  u h ð  x ð n ; h Þ ; t  Þ ; ~ e q ð n ; h ; t  Þ ¼  q ð  x ð n ; h Þ ; t  Þ   q h ð  x ð n ; h Þ ; t  Þ ; ~^ e u ð 1 ; h ; t  Þ ¼  ^ e u ð  x ð 1 ; h Þ ; t  Þ ;  ~^ e q ð 1 ; h ; t  Þ ¼  ^ e q ð  x ð 1 ; h Þ ; t  Þ ; V  ð n Þ ¼  v  ð  x ð n ; h ÞÞ ;  W  ð n Þ ¼  w ð  x ð n ; h ÞÞ ; we obtain the LDG orthogonality condition (2.13) on the referenceelement [  1,1] h 2 Z   1  1 ð ~ e u Þ tt  Vd n þ  c  2 Z   1  1 ~ e q V  0 d n    c  2 ~^ e q ð 1 ; h ; t  Þ V   ð 1 Þþ  c  2 ~^ e q ð 1 ; h ; t  Þ V  þ ð 1 Þ ¼  0 ;  ð 2 : 15a Þ h 2 Z   1  1 ~ e q Wd n þ Z   1  1 ~ e u W  0 d n   ~^ e u ð 1 ; h ; t  Þ W   ð 1 Þþ ~^ e u ð 1 ; h ; t  Þ W  þ ð 1 Þ ¼  0 :  ð 2 : 15b Þ If the exact solutions  u  and  q  are analytic functions, then we can ex-pand the local errors as a Maclaurin series with respect to  h  as e u ð  x ð n ; h Þ ; t  Þ ¼ X 1 k ¼ 0 U  k ð n ; t  Þ h k ; e q ð  x ð n ; h Þ ; t  Þ ¼ X 1 k ¼ 0 Q  k ð n ; t  Þ h k ;  ð 2 : 16 Þ where  U  k (  , t  ) and  Q  k ð ; t  Þ 2 V  kh  are obtained by applying the chainrule as U  k ð n ; t  Þ ¼  1 k ! d k e u dh k  ð  x ð n ; 0 Þ ; t  Þ ¼  1 k ! X ki ¼ 0 n i 2 i ki  @  i x @  k  ih  e u ð 0 ; t  Þ ;  ð 2 : 17 Þ Q  k ð n ; t  Þ ¼  1 k ! d k e q dh k  ð  x ð n ; 0 Þ ; t  Þ ¼  1 k ! X ki ¼ 0 n i 2 i ki  @  i x @  k  ih  e q ð 0 ; t  Þ ;  ð 2 : 18 Þ where the binomial coefficient  ki    is defined by  ki   ¼  k ! i ! ð k  i Þ !  for 0 6 i 6 k . In the remainder of this paper we will omit the  unless we feelit is needed for clarity. 132  M. Baccouch/Comput. Methods Appl. Mech. Engrg. 209–212 (2012) 129–143  3. Local error analysis Here we investigate the superconvergence properties for theLDG method applied to the wave Eq. (2.1). In particular, weshow that the LDG solution is  O ( h  p +2 ) superconvergent at the(  p  + 1)-degree right-Radau polynomial and the solution’sderivative is  O ( h  p +2 ) superconvergent at the (  p  + 1)-degree left-Radau polynomial. The local superconvergence results are provedand the global superconvergence results are confirmednumerically.For simplicity, we present an error analysis on the element I   = [0, h ]. We consider the problem (2.1a) on [0, h ] subject to the ini-tial condition (2.1b), where  x 2 [0, h ] and the Dirichlet boundaryconditions  u (0, t  ) =  u l ( t  ),  u ( h , t  ) =  u r  ( t  ),  t  2 [0, T  ]. The proof of our re-sults with the mixed boundary conditions is the same and isomitted.Since we consideron element,the fluxes (2.6g)usingthe Dirich-let boundary conditions with  d  ¼  ph  become ^ u h ð 0 ; t  Þ ¼  u l ð t  Þ ;  ^ u h ð h ; t  Þ ¼  u r  ð t  Þ ;  ^ q h ð 0 ; t  Þ ¼  q h ð 0 ; t  Þ ; ^ q h ð h ; t  Þ ¼  q h ð h ; t  Þ þ  ph ð u h ð h ; t  Þ  u r  ð t  ÞÞ ; where  u h (  x , t  ) and  q h (  x , t  ) are the LDG solutions on [0, h ]. Using themapping of [0, h ] onto [  1,1] given by (2.14), we have ^ e u ð 1 ; h ; t  Þ ¼  ^ e u ð 1 ; h ; t  Þ ¼  0 ;  ^ e q ð 1 ; h ; t  Þ ¼  e q ð 1 ; h ; t  Þ ; ^ e q ð 1 ; h ; t  Þ ¼  q ð 1 ; h ; t  Þ  ^ q h ð 1 ; h ; t  Þ¼  q ð 1 ; h ; t  Þ  q h ð 1 ; h ; t  Þ   ph ð u h ð 1 ; h ; t  Þ  u 1 ð t  ÞÞ¼  e q ð 1 ; h ; t  Þþ  phe u ð 1 ; h ; t  Þ : The LDG orthogonality conditions (2.15) for the local errors  e u  and e q  on the element [0, h ] can be simplified to h 2 Z   1  1 ð e u Þ tt  Vd n þ  c  2 Z   1  1 e q V  0 d n   c  2 e q ð 1 ; h ; t  Þ þ  phe u ð 1 ; h ; t  Þ  V   ð 1 Þþ  c  2 e q ð 1 ; h ; t  Þ V   ð 1 Þ ¼  0 ;  ð 3 : 1a Þ h 2 Z   1  1 e q Wd n þ Z   1  1 e u W  0 d n  ¼  0 :  ð 3 : 1b Þ Since we consider only one element, we will omit the  ± , for instance, V  + (  1) =  V  (  1) and  V   (1) =  V  (1). In the following theorem, we state and prove the main result of this section. It gives the superconvergence properties for the LDGmethod applied to the wave Eq. (2.1a), subject to the initial condi-tions (2.1b) and to the Dirichlet boundary conditions (2.1c).  Theorem 3.1.  Let (u,q) and (u h, q h ), respectively, be the solutions of  (2.3)  and  (2.6)  in I =[0,h] subject to ( 2.1b ) and ( 2.1c ). Then, the local finite element errors can be written as −1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1−1−0.500.511.52−1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1−0.015−0.01−0.00500.0050.010.0150.020.0250.03−1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1−2−101234x 10 −4 −1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1−1.5−1−0.500.511.522.533.5x 10 −6 Fig. 1.  Zero-level curves of   e u (  x ,5) for Example 5.1 using  P   p  with  p  = 1  4 on uniform meshes having  N   = 12 elements (upper left to lower right). M. Baccouch/Comput. Methods Appl. Mech. Engrg. 209–212 (2012) 129–143  133
Related Documents
View more...
We Need Your Support
Thank you for visiting our website and your interest in our free products and services. We are nonprofit website to share and download documents. To the running of this website, we need your help to support us.

Thanks to everyone for your continued support.

No, Thanks