Convergence of discontinuous Galerkin schemes for front propagation with obstacles
We study semi-Lagrangian discontinuous Galerkin (SLDG) and Runge-Kutta discontinuous Galerkin (RKDG) schemes for some front propagation problems in the presence of an obstacle term, modeled by a nonlinear Hamilton-Jacobi equation of the form , in one space dimension. New convergence results and error bounds are obtained for Lipschitz regular data. These “low regularity” assumptions are the natural ones for the solutions of the studied equations. Numerical tests are given to illustrate the behavior of our schemes.
Key words and phrases:Hamilton-Jacobi-Bellman equations; discontinuous Galerkin methods; level sets; front propagation; obstacle problems; dynamic programming principle; convergence
In this paper, we establish convergence of a class of discontinuous Galerkin (DG) methods for the one-dimensional Hamilton-Jacobi (HJ) equation below, hereafter also called the “obstacle” equation,
with periodic boundary conditions on and a constant . In (1), the function is called the “obstacle” function.
It is well known that taking constraints in optimal control problems is not an obvious task. Within the viscosity theory, it is possible to devise schemes for obstacle equations such as (1). However a monotonicity condition is needed in general for proving the convergence of the scheme. The monotonicity condition can yield a convergence proof of one-half order in the mesh size  (see also  for more specific partial differential equations (PDEs) with an obstacle term and related error estimates). However, a serious restriction of such monotonicity condition is that the schemes become at most first order accurate for smooth solutions or in smooth regions, thus making the schemes highly inefficient for practical computation. On the other hand, it is very difficult to show convergence of formally higher order schemes, such as the ones studied in this paper, when the solution is not regular enough.
In our previous work , we proposed a class of Runge-Kutta DG (RKDG) methods adapted to front propagation problems with obstacles. The DG methods under consideration were originally devised to solve conservation laws, see for example the review paper . As for the DG-HJ solvers, in [19, 23], the first efforts relied on solving the conservation law system satisfied by the derivative of the solution. See also  for an adaptive version of this scheme. In , a DG method for directly solving the Hamilton-Jacobi equation was developed and was later generalized to solve front propagation problems  and obstacle problems . Other direct DG solvers include the central DG scheme  and the local DG scheme . The schemes proposed in  feature a simple treatment of the obstacle functions. Stability analysis is performed with forward Euler, a Heun scheme and a TVD Runge-Kutta third order (TVD-RK3) time discretization using the techniques developed in Zhang and Shu .
On the other hand, the semi-Lagrangian DG (SLDG) methods were proposed in [27, 28, 26] to compute incompressible flow and Vlasov equations, as well as in  for some general linear first and second order PDEs. The advantage of SLDG is its ability to take large time steps without a CFL restriction. However, it is difficult to design SLDG methods for nonlinear problems (some SL schemes for Hamilton-Jacobi-Bellman equations were proposed in , but without convergence proof). For general SL methods, we refer to the works of Falcone and Ferretti [14, 15, 16], see also the textbook . There is also a vast literature on high order finite difference schemes for solving HJ equations, see, e.g. [25, 1, 21, 22].
Beyond the scope of the present paper, yet of great interest, we also mention the second order PDE with obstacle terms, such as , with . This is the case of the so-called “American options” in mathematical finance . Explicit schemes were proposed and proved to converge within the viscosity theory (see for instance ) yet with a reduced rate of convergence. Variational methods for nonlinear obstacle equations can also be devised  but will lead in general to nonlinear implicit schemes which are computationally more demanding.
The scope of the present paper is to study convergence of the SLDG and RKDG schemes for the obstacle problem (1). The main challenges include the low regularity of the solution and the nonlinear treatment needed to obtain the obstacle solution. Due to the fully discrete nature of the method, traditional techniques for obtaining semi-discrete error estimates of DG methods for hyperbolic problems cannot directly apply here. Therefore, fully discrete analysis is necessary. Fully discrete analysis of RKDG methods for conservation laws has been performed in the literature. In , error estimates for RKDG methods for scalar conservation laws were provided for smooth solutions. In [11, 31], discontinuous solutions with RK2 and linear polynomials and RK3 time discretizations with general polynomials were studied for linear conservation laws. However, to our best knowledge, convergence results for second and higher order DG schemes solving nonlinear hyperbolic equations with irregular solutions are not available.
As a consequence of our results, assuming that is a space step and a time step, we shall show error bounds of the order of for the SLDG schemes (under time stepping , larger time step can be taken with a lower convergence rate), and of order for the RKDG schemes (under time stepping ). We will need a natural “no shattering” regularity assumption on the exact solution that will be made precise in the sequel, otherwise we will typically assume the exact solution to be Lipschitz regular and piecewise regular for some for SLDG (resp. for RKDG).
The main idea of our proof is based on the dynamic programming principles as illustrated below. The viscosity solution of (1) also corresponds to the following optimal control problem:
The function is also the solution of the Bellman’s dynamic programming principle (DPP): for any ,
If we denote
then the DPP implies in particular for any and :
Using formula (2), we can see that when and are Lipschitz regular, then is also Lipschitz regular in space and time, and in general no more regularity can be assumed (the maximum of two regular functions is in general no more than Lipschitz regular). When is a discontinuous function (otherwise piecewise regular), formula (2) implies also some regularity on the solution .
The rest of the paper is organized as follows. In Section 2, we describe the DG methods for the obstacle equations. In Section 3, we collect some lemmas that will be used in our convergence proofs. Section 4 and Section 5 are devoted to the convergence analysis of SLDG and RKDG methods, respectively. In both sections, we will first establish error estimates for SLDG and RKDG schemes for linear transport equations without obstacles. Then we will use DPP illustrated above to prove convergence of the numerical solution in the presence of obstacles. Numerical examples are given in Section 6. We conclude with a few remarks for the multi-dimensional case in Section 7.
2. DG schemes for the obstacle equation
In this section, we will introduce the SLDG and RKDG methods for the obstacle equation (1). For simplicity of discussion, in the rest of the paper, we will assume to be a positive constant.
Let be a set of intervals forming a partition of . We denote where is the length of the interval . For a given integer , let be the DG space of piecewise polynomial of degree at most on each interval :
To introduce the methods for the obstacle problems, we follow two steps. Firstly, we describe the DG solvers for the linear advection equation .
To advance the numerical solution in one step from to , we consider one of the two DG methods described below.
where is the projection onto the space . We shall denote this SLDG solver by .
RKDG Scheme (by TVD-RK3 time stepping): find , , , such that
We shall denote this RKDG solver by .
After introducing DG schemes for the linear transport problem, for the obstacle equation (1), we shall consider two approaches: one by using the projection
where or and
The idea of the formulation above is to try to follow the relation (4), and the projection step is to project the function into the piecewise polynomial space .
Unfortunately, the scheme (8) is difficult to implement, because we need to compute the maximum of two functions, which requires locating the roots of . Another more practical approach is to define as the unique polynomial in such that
where are the Gauss-Legendre quadrature points on the interval , and are the corresponding quadrature weights. Those schemes were studied in details and stability was established in . We shall see in later sections that definitions (8) or (10) lead to similar error estimates, although the second approach is much easier to implement.
Finally, we remark that in , forward Euler and TVD-RK2 temporal discretizations are also considered. However, the stability restriction for the time step for the forward Euler method is rather severe as , and the stability proof of TVD-RK2 only works for piecewise linear polynomials. Therefore, in this paper, we will only consider TVD-RK3 time discretizations.
In this section, we collect some lemmas which will be used in our convergence proof, and discuss properties of the obstacle solutions. Here and below, we use (possibly with subscripts) to denote a positive constant depending solely on the exact solution, which may have a different value in each occurrence.
Let us introduce, for , the following function sets:
where denotes the -th derivative almost everywhere, and
The following pseudo-norm definition will also be used:
In particular, using the Gauss-Legendre quadrature rule, for any we have . From this point on, we will use to denote , and to denote for a given domain D.
3.1. Properties of projections and the obstacle function
For the RKDG method, it is necessary to consider the following Legendre-Gauss-Radau projection . For any function , , and for any element , it holds that
In the lemma below, we will first establish the projection properties for functions in the space .
Let and let be in : is Lipschitz continuous, piecewise , with at most non regular points. Then there exists a constant , depending only on such that
where the projection can be either or , and is defined as
Using the property of the projections , on each regular cell , we have
since on . Therefore,
On the other hand, on the intervals where is not regular, we have
because is Lipschitz continuous. Hence summing up on the “bad” intervals (with at most such intervals),
Summing up the bounds with bad intervals and good ones, we prove the desired result. ∎
We now state a similar estimate for the norm:
Assume that belongs to , , then we have
where or , is defined as in (12), and the constant depends only on and .
On each regular cell , we have ,
Then using the fact that , we obtain that
On the other hand, when the interval is such that contains a non-regular point, we can write
This concludes our proof. ∎
We now turn to some estimates related to the obstacle function .
Assume that for some integer . Let . There exists a constant such that
Denote the set of local maximum points of , if , then we see that . Furthermore
since there are at most local maxima.
Consider the case when contains at least a local maxima of . Assuming that is small enough we can assume that is the only local maximum of on the interval , so that . Then,
since is twice differentiable at .
For the second estimate (15), we make use of a minimal covering of the set
using mesh intervals. The length of this covering is bounded by for each maximum point, since in order to cover any interval we may need two more mesh intervals , of length than the minimum required length . Overall the length of the total covering is bounded by , hence we obtain the desired result.
3.2. Properties of the obstacle solutions
We shall impose some restrictions on the regularity of the obstacle solutions as described below.
[“no shattering” property] For a given , we will say that the exact solution of the problem (1) is “not shattering” if there exists some and constants such that the exact solution satisfies, for any , .
Recall that the exact solution satisfies . Therefore if is not shattering, it implies that has bounded number of zeros (since otherwise would have an unbounded number of singularities).
A typical example where shattering occurs (therefore not satisfying definition 3.1), can be constructed as follows. Suppose the domain contains the interval , let and on , together with velocity constant . The function is of class on the interval . Notice then, for and , and the exact solution is therefore for and . So at time , (for ). This function has an infinite number of non regular points in the interval , and therefore does not satisfy the “no shattering” property at time .
It is not easy to state precise conditions on the initial data and to ensure that the no shattering property will be satisfied. Mainly, , should have only a finitely bounded number of zeros, as is detailed below. However, it is clear that shattering will not occur for generic data and . This definition still allows for a finite (bounded) number of singularities in , as is generally the case when taking the maximum of two regular functions. Finally we also give an example (see Lemma 3.5 below) where we can prove that the “no shattering” condition is satisfied.
Assume that, for a given ,
|and that, ,|
the bound being independent of .
Then there exists constants such that, , belongs to (with , , and that are independent of .)
First, for , if we denote (resp. ) to be the number of local maxima (resp. minima) of , then we have . This is because the Lipschitz constant of is bounded by by an elementary verification. Each local maxima of may develop into two singularities in , each local minima may develop into one singularity in , and each singular point of may continue to be a singularity in . Hence the total number of non-regular points in will be bounded by . The bound of the derivative can also be obtained easily.
Because the exact solution is given by , the Lipschitz constant of is therefore bounded by .
On the other hand, the number of singular points of is bounded by the sum of the number of singular points of the function , plus the number of zeros of the same function, which is assumed to be bounded independently of . Hence the the number of singular points of is bounded independently of .
Finally the bound on the -partial derivative , in the regular region of , is easily obtained, because is then locally one of the two functions or . ∎
4. Convergence of the SLDG schemes
In this section, we will provide the convergence proof of the SLDG scheme. In particular, we will proceed in three steps. First, we will establish error estimates of the SLDG methods for the linear transport equation
4.1. Convergence of the SLDG scheme for the linear advection equation
We first consider the linear equation , for which
We denote , and we define the numerical solution of the SLDG method at to be . In particular, the scheme writes: initialize with , and for .
We consider , . If , then we have for all such that ,
for some constant independent of , and is defined in (12).
If , then , and due to Lemma 3.1, we have:
for some constant independent of . Hence
where we have used the fact in the fourth row, and the periodic boundary conditions in the last row. As for the initial condition, we have
Finally we obtain for any given the existence of a constant (independent of and of ), such that
for all such that . ∎
Therefore, if , then , and we need for the convergence. If otherwise , it holds and we would only need for the convergence.
4.2. Convergence of the first SLDG scheme in the obstacle case
Assume that the exact solution is not shattering in the sense of Definition 3.1, for some integer . Then the following error bound holds:
for some constant independent of . In particular,
If , then :
If , and if , then :
Defining as in , and assuming , the error is bounded by . Therefore the optimal estimate is obtained when , or , and the error is of order .
4.3. Convergence of the SLDG scheme defined with Gauss-Legendre quadrature points
Let and assume that the exact solution is not shattering in the sense of Definition 3.1. The following error bound holds:
for some constant independent of . In particular,
If , then :
If and , then :