Superconvergence of the local discontinuous Galerkin method for nonlinear convection-diffusion problems

In this paper, we discuss the superconvergence of the local discontinuous Galerkin methods for nonlinear convection-diffusion equations. We prove that the numerical solution is \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$(k+3/2)$\end{document}(k+3/2)th-order superconvergent to a particular projection of the exact solution, when the upwind flux and the alternating fluxes are used. The proof is valid for arbitrary nonuniform regular meshes and for piecewise polynomials of degree k (\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$k\geq1$\end{document}k≥1). The numerical experiments reveal that the property of superconvergence actually holds true for general fluxes.


Introduction
In this paper, we discuss the nonlinear convection-diffusion equations given by with the periodic boundary condition, where ν >  is a constant. We study the superconvergence of the local dicontinuous Galerkin (LDG) solutions towards a particular projection of the exact solution.
The high-order numerical methods have been applied in a variety of fields [-]. The LDG method is one of those numerical methods, which were first constructed by Cockburn and Shu and motivated by Bassi and Rebay [, ] to solve the convection-diffusion equations. Since then, the LDG method has been used to solve the time-dependent equations with high spatial derivatives, such as the Korteweg-de Vries (KdV) equations [], time-dependent fourth-order problems [] and the general fifth-order KdV equations []. See more details in []. We now state some theoretical results, which represent the crucial technique to treat the nonlinear parts of the equations. In [], Zhang and Shu study the error estimate of the discontinuous Galerkin (DG) method with second-order Runge-Kutta time discretization. They obtain the optimal error estimate of the (k + )th order for upwind numerical fluxes and a suboptimal error estimate of the (k + /)th order for general monotone fluxes, where k is the order of the piecewise polynomial space. The proof holds true for arbitrary meshes under the reasonable assumptions. Then Zhang and Shu extend the results in [] to the third-order TVD Runge-Kutta time discretization case, which is more popular in the computation []. In [], Wang and Shu obtain the optimal error estimate of the LDG method with implicit-explicit time-marching for nonlinear convection-diffusion problems, when the fluxes are chosen to be the general monotone fluxes and alternating fluxes. In the above three papers, the Taylor expansion and an a priori assumption are used to estimate the nonlinear parts.
We would like to mention the superconvergence results for DG and LDG methods. In [], Cheng and Shu study the superconvergence of the (k + /)th order of the DG solution towards a particular projection for linear conservation laws. The limitations of [] are that the proof is only valid for uniform meshes and linear piecewise polynomial space. Cheng and Shu overcome this limitation in [], which implies that the result in [] holds true for arbitrary meshes and kth-order finite element spaces. Cheng and Shu also extend the result to the linear convection-diffusion problems. For the linear equations with highorder spatial derivatives, Hufford obtains the superconvergence of the (k + /)th order for linear KdV equations [] and Meng gets the same result for the linear fourth-order problems []. But for above linear problems, the numerical experiments imply that the numerical solution is superconvergent to the exact solution at a rate of the (k + )th order. It is highly nontrivial to obtain this half-order increase theoretically. For linear conservation laws and linear parabolic equations, Yang and Shu use a new technique to carry out the optimal order of the superconvergence [, ]. In addition, they prove that DG and LDG solutions are (k + )th-order superconvergent to the exact solutions at Radau points. In [-], Cao and Zhang present another framework to demonstrate the superconvergence at Radau points for linear -D and -D hyperbolic problems and -D linear parabolic problems. The first superconvergence proof with the DG method for nonlinear conservation laws is obtained in [], when the upwind fluxes are used, under the condition that the absolute value of the convection term f has a positive low bound. In this paper, we obtain a similar result for the nonlinear convection-diffusion problems, when the upwind fluxes and alternating fluxes are used under the assumption that |f | ≥ . Due to the character of the LDG method, there is no need of a strict positive bound of the absolute value of the convection term. In [], Cao obtains the superconvergence of DG methods based on upwind-biased fluxes for -D linear hyperbolic equations. Guo and Yang show the DG solution is (k + )th-order accurate at the downwind points and (k + )th-order accurate at all the other downwind-biased Radau points in [].
The outline of this paper is as follows. In Section , we present the semi-discrete LDG schemes for nonlinear convection-diffusion problems. In Section , we state the main proofs of our theorems. Some numerical experiments are presented in Section , and in Section , we give the conclusion and our future work. Finally, we give a proof of a lemma in the Appendix.

The local discontinuous Galerkin method for nonlinear convection-diffusion equations 2.1 The local discontinuous Galerkin method
In this subsection, we will present the semi-discrete LDG method for equation (.). We first divide the computational domain = [, π] into N subintervals. We have We denote each subinterval by I j = [x j-/ , x j+/ ] and the union of all I j by I h . Let h j = x j+/x j-/ be the length of the subinterval and h = max ≤j≤N h j . We denote the left and right limits of the function v h at the element interface x j+/ by (v h ) -j+/ and (v h ) + j+/ , respectively. We also set j+/ and denote the center of I j by x j = (x j+/ + x j-/ )/. We would like to assume our mesh is regular, which means that there exists a constant λ >  such that λh < h j .
We choose the finite element space as the kth-order piecewise polynomial space that is denoted by where P k (I j ) is the space of polynomials of degree at most k on I j . In addition, we define the broken Sobolev space on = [, π] as Before we construct the LDG method, we need an auxiliary variable q, so we rewrite equation (.) as a first-order system. We have Then the semi-discrete LDG scheme is formulated as follows: find u h , q h ∈ V k h such that for any wheref (a, b) is usually a monotone flux, which satisfies: ) is a Lipschitz continuous function in both arguments.
•f (a, b) is a nondecreasing function in a and a nonincreasing function in b. Hereû h ,q h are the alternating fluxes, of which we have two choices: A new technique is necessary when the sign of the derivative of the convection term changes. Hence we only consider the case of a sign-preserving derivative, which implies that we could use the upwind flux. We havê

Norms
In this subsection, we give some norms used in this paper. We denote the standard L  norm on I j by · I j . Then the norm of Sobolev space H l (I j ) is defined as where l is a natural number and D α is the αth-order spatial derivative operator. We denote the norm of H l h by For convenience, we set · = ·  . Moreover, the L ∞ norm of the whole computational domain is where u ∞,I j is the standard L ∞ norm on I j . The norm on the boundary of I j is defined as Then we have

Properties of the finite element space
We first present the Gauss-Radau projections, from H k+ h into V k h , which are denoted by In addition, we define P + h u as for any I j .
We denote the projection error (u -Ph u) or (u -P + h u) by η u . Thanks to the standard approximation theory [], it is easy to obtain the following approximation property. Suppose u(x) ∈ H k+ (I j ) and is sufficiently smooth. Then we have where C  is a positive constant independent of h j .
Here and below, we use C  to denote the corresponding constant of the estimate for projection errors. Now we turn to the three inverse properties of the finite element space there exists a positive constant μ independent of h and j, such that For more details on these inverse properties, we refer the reader to [].

Initial projection
In order to complete the proof, an initial condition compatible with the superconvergence order would be constructed with care. Fortunately, a (k + /)th-order initial condition P h presented in [] is valid in our proof, that is, for any function u, for any I j . Also, we require We would like to remark that P h exists and is unique. Moreover, we have the following estimate:

Superconvergence
In this section, we will give the main proof of superconvergence. Similar to [], we assume that the exact solution u(x, t) is sufficiently smooth. We have where S is a constant independent of t and h. Also, the flux function f is smooth enough, for example f ∈ C  , and its derivatives are bounded on R. We have We would like to remark that this assumption is reasonable with the original or a suitably modified flux f , if we only consider the smooth solution. See for more details [].
Without loss of generality, we assume f >  and ν = . Then the fluxes are chosen aŝ For ease of notation, we denote the error between the exact solution u and the numerical solution u h by e u = uu h . Also, we set Next, we will introduce three operators, namely, When the fluxes (.) are used, we rewrite (.a)-(.b) as For u, q are exact solutions, we get a similar system to (.a)-(.b), which is Then we obtain the error equations Due to the properties of Gauss-Radau projections and summing (.a)-(.b) over j, we obtain where K = N j= K j and H ± = N j= H ± j . Now we turn to an investigation of the properties of H ± (·, ·) and K(·, ·). By integrating by parts, we get Under the periodic conditions, we have the following equation: whose proof is straightforward. For the nonlinear part, the derivation of the estimate is a little involving. We would like to remark that the following estimate inequality is not the final form, since we will adjust it further to obtain the suitable form in different proofs.
Lemma  Under the assumptions (.) and (.), we have the following estimate of the nonlinear part. For any v h ∈ V k h , where C f is a constant dependent on |f |, |f |, μ and the exact solution u, but independent of h.
Proof We begin with using the second-order Taylor expansion with respect to the variable u. We have where We remark that we have removed the subscript (j + /) in (.) for notational convenience.
Noting that (ηu ) j+/ = , we divide (K(u, v h ) -K(u h , v h )) into three parts. We have We will estimate these three parts below: • The estimate of |  |. Due to property (.a)-(.b), we have For |f | is bounded, it is easy to show that |f (u)f (u j )| ≤ C f h. Then we employ the Cauchy-Schwarz inequality to obtain Proceeding as in the estimate of |  |, we split the integration into two parts. We have Noting that with straightforward application of the Cauchy-Schwarz inequality and using properties (.a) and (.c), we obtain Applying the Cauchy-Schwarz inequality, we obtain In accordance with properties (.b) and (.c), we have Collecting (.), (.) and (.), we obtain the estimate (.).
We would like to remark that the estimate is also valid on each element I j , though we only present the result on the whole computational domain . Now we move on to the error estimate of semi-discrete LDG schemes.
Lemma  Suppose that u, q are the exact solutions of system (.a)-(.b), which are sufficiently smooth, and that u h , q h are the solutions of (.a)-(.b). Under assumptions (.) and (.), we have the following estimate: Proof Taking (v h , w h ) = (ξ u , ξ q ) and adding up (.a) and (.b), we have We deduce from property (.) that On the other hand, plugging v h = ξ u back into (.), applying the Cauchy-Schwarz inequality to the integration term and using property (.a), we obtain Now we turn to the estimate of the last term on the right hand side of (.), as follows: In the first step, we use the triangle inequality and property (.). The proof of the last step is an application of Young's inequality. Then we employ property (.) and Young's inequality to get Finally, applying property (.), the Cauchy-Schwarz inequality and Young's inequality, we have Collecting (.) and (.), we arrive at (.). We thus finish our demonstration.
To estimate e u and e q , we would like to follow [] to make the a priori assumption that, if h is sufficiently small, we have It follows from property (.c), property (.) and the triangle inequality that We will verify (.) at the end of this paper.
If assumption (.) is true, then there exists a positive constant δ such that By using the Gronwall inequality and the triangle inequality, we obtain the error estimate as follows: where C is a constant associated with |f |, |f |, μ, the exact solution u and the final time T.
In addition, we give a lemma to estimate (e u ) t , which will be proved in the Appendix.
Lemma  Under the assumptions (.), (.) and (.), we have the following estimate: Next we will present a crucial lemma which will be used to derive the superconvergence.
Lemma  Suppose assumptions (.), (.) and (.) are true. Then we have the following inequalities: where C s is a constant independent of h.
Proof The first inequality is the same as Lemma . given in [], so only the second inequality will be proved. Applying property (.) to equation (.a), we get According to properties (.) and (.), we have Using the Cauchy-Schwarz inequality and property (.) yields where L(x) is the kth-order Legendre polynomial on [-, ]. Then we obtain We are now in a position to prove our theorem.
Theorem  Suppose that u, q are the exact solutions of (.a)-(.b), which are sufficiently smooth, and that u h , q h are the solutions of (.a)-(.b). We also assume that f ∈ C  and |f |, |f |, |f | are bounded on R.
where the positive constant C * is independent of h, but maybe depends on u, f and T.
Proof Recall inequality (.). To obtain the half-order increase, we shall use Lemma  to estimate the first term on the right side of (.). We have (.) Setξ j u =  h j I j ξ u dx. Then, according to the orthogonality of the Gauss-Radau projection, we get It is a simple application of the Cauchy-Schwarz inequality to obtain the second inequality. The third one follows from the Poincaré inequality. Applying Lemma , the Cauchy-Schwarz inequality and property (.), we arrive at the last two inequalities.
Collecting (.), (.) and (.) and using Young's inequality, we obtain The estimates of the remaining two arguments are quite similar to (.), so we only present the results as follows: According to Lemma  and estimate (.), we have Recalling estimate (.), we obtain where C,C are constants independent of h.
Integrating with respect to t, it follows from Gronwall's inequality and the estimate of the initial projection that Finally, we will verify the a priori assumption (.) to complete our demonstration. We first mention that there exists a positive h  , for any h < h  , such that C * h k+/ <   h  and C IP h k+/ <   h  , where C is the constant in (.) and C IP is the constant in (.). Then, when t = , for any h < h  , we have We now define For M not empty, we denote the supremum value of M by t sup . If t sup < T, it follows from the continuity of ξ u (·, t) that Then we actually have which is a contradiction to (.). Hence, we have t sup = T, which justifies the a priori assumption (.). Thus we finish our proof.

Numerical experiments
In this section, we will give some numerical results to support our theorems. In all experiments, the time discretization is the third-order IM-EX Runge-Kutta scheme

Example 1
We first consider the following equation with periodic boundary condition: u(x, ) = sin(x).

(.)
In this case, f = u  > , which implies that we can use upwind fluxes. In Table  and Table , we present the L  errors of e u and ξ u and their orders on a nonuniform mesh, which is a % random perturbation of the uniform mesh, at the final time T =  in the P  piecewise polynomial case and the final time T = . in the P  piecewise polynomial case, respectively.

Example 2
We take the equation of which the flux function changes its sign on the computational domain. Hence, we use the Godunov flux in this example. The numerical results on the nonuniform mesh, which is a % random perturbation of the uniform mesh, are presented by Table  and Table , which imply that the superconvergence property is still valid in the case that the flux function is not sign preserving.

Example 3
In this example, we take an equation with a non-polynomial flux function. We have The boundary condition is a periodic boundary condition. The mesh is also a % random perturbation of the uniform mesh. It results from Table  and Table  that the superconvergence property is true for a strong nonlinear flux function.

Conclusion
In this paper, we investigate the superconvergence of the LDG method for nonlinear convection-diffusion problems. The order of the superconvergence of the LDG method with P k (k ≥ ) piecewise polynomial as the finite element space is proved to be the (k + /)th-order when the fluxes are upwind fluxes and alternating fluxes. The numerical experiments demonstrate that the superconvergence property is valid for general flux functions.  Future work includes the study of superconvergence of the LDG method for the nonlinear equations with high-order spatial derivatives in -D. The superconvergence properties of general monotone numerical flux will also be considered. where a)-(A.b) and applying property (.) gives We now turn to an estimate of |N L(u, u h , (ξ u ) t )|. The process is similar to Lemma , except for the fact that we use the second-order Taylor expansion, whose remainder is of the integral form, to deal with ∂ t (f (u)f (u h )). We have We are now ready to estimate each part: • Noting the properties of the Gauss-Radau projection, we have • If we do what we did in (A.) and (A.), we obtain • Dividing the integration into two parts gives In the second inequality, we use the Cauchy-Schwarz inequality. Using properties (.a) and (.b), we get the final step. • Similar to the process of (A.), we obtain Following estimate (.), Gronwall's inequality and the triangle inequality, we complete the proof.