Original articlesImproving the accuracy of LDG approximations on coarse meshes
Introduction
Finding a mesh that produces an accurate approximation is an old problem in finite element analysis[[28], [32], [7], [17], [18]]. Most of the work in mesh optimization is guided by the principle of discrete potential energy minimization. The optimization problem assumes the energy to be a function of the nodes coordinates and seeks the optimal node position among all the possible triangulations of the computational domain , see the pioneering work of Carroll and Bakker [7], McNeice and Marcal [28], and Felippa [[17], [18]]. However this is an extremely difficult problem since it must consider a variable number of nodes as well as a variable mesh connectivity. In practice, the problem is simplified by restricting the search in the set, , of all the triangulations of having the same connectivity, hence the same number of nodes. This procedure is referred by some authors as -adaptivity.
In [28], McNeice and Marcal derived some necessary conditions to attain the minimum value of the discrete energy and provided a rudimentary algorithm, suitable for meshes with a small set of nodes. A close related work was also presented by Carroll and Bakker [7]. In [17], Felippa suggested, for the first time, the incorporation of more sophisticated unconstrained minimization techniques for an automatic grid optimization process. However, the application of efficient optimization techniques most often requires the explicit knowledge of the gradient’s minimization functional or its approximation. Thus, in grid optimization most of the theoretical research has focused on deriving an explicit formulation of the energy gradient. For instance, Delfour, Payre and Zolésio [15] presented an explicit expression for the directional error gradient for the continuous finite element method; Oliveira [30] showed the existence of optimal meshes and carried out a theoretical analysis of the differentiability of the error with respect to the node coordinates in [29]. More recently Tourigny and Hülsemann [36] analyzed a global minimization procedure based on local Dirichlet variational problems; and derived an explicit expression for the gradient of the local problems.
Most of the research in -adaptive methods has been focused on time dependent problems. In that context, the technique is better known as moving finite element method. Simply explained, the original (or physical) problem is coupled with an additional PDE controlling the position of the grid points over time; and the approximation on the new grid is obtained by interpolation or by performing a recalculation using the updated mesh. Central to this technique is the appropriate definition of a monitor function which couples the mesh and the physical solution. For instance, Huang and Rusell [20] proposed a monitor function based on the gradient of the solution for two dimensional problems. Since the rationale is to concentrate grid points in regions where the solution has steep gradients; monitor functions based on a posteriori error indicators have also been considered by Cao, Huang and Russell [6]; and by Lang et al.[23], in combination with a local adaptive strategy. Another interesting approach, was proposed by Li and collaborators [[26], [24]] as an extension of Dvinsky’s work [16] on grid generation using the harmonic function theory. For a comprehensive literature review on the moving mesh technique the reader is referred to Huang [19] and Tang [34]; and references therein.
Discontinuous Galerkin (dG) method is, by now, a well known spatial discretization technique that has been, extensively, used in the approximation of partial differential equations. The Local Discontinuous Galerkin (LDG) method, originally introduced by Cockburn and Shu [13], is one among several dG methods, that has been, successfully, applied to a wide variety problems. The convergence and a priori error analysis for the LDG method applied to elliptic problems can be found in [[10], [31], [12]]. Due to its excellent stability and approximation properties, see [8] for a comparative study of several dG methods in 2D; this work concentrates on the LDG method. However, we must point out that the proposed technique could also be applied to other dG methods.
To put the present work into perspective, let us briefly recall the scarce work on -adaptive techniques applied to dG discretizations. In [25] Li and Tang extended their previous work on moving meshes to nonlinear conservation laws using the Runge Kutta Discontinuous Galerkin (RKDG) method, [14], as time–space discretization. They considered a monitor function based on a simple a posteriori error estimator and a suitable interpolation operator. Using the framework developed in [20], Mackenzie and Nicola [27] proposed another technique for the Hamilton–Jacobi equations and discontinuous spatial discretizations which preserves the conservation properties of the problem. Recently, using the Interior Penalty (IP) method [2] as spatial discretization, Uzunca, Karasözen and Küçükseyhan [37] presented a moving mesh technique based on the rezoning approach [21] to solve nonlinear PDEs with traveling waves. Antonietti and Houston [1], proposed an iterative moving mesh technique for steady state convection–diffusion problems where the monitor function is based on goal-oriented a posteriori error indicators. Their numerical experiments showed that more accurate results could be obtained when an -adaptive strategy is applied after improving the initial coarse mesh using their moving mesh algorithm.
Typically, an accurate finite element approximation is obtained after performing few steps of an -adaptive technique; hence increasing the total number of cells. This can be a major drawback for the LDG method who suffers from having a larger stencil compared to other dG methods having a compact stencil.
The aim of this work is to obtain a more accurate approximation for a given mesh connectivity; by finding an optimal mesh (or closed to it) that minimizes the discrete potential energy functional of the LDG method. Therefore, since the mesh connectivity remains fixed no increase in memory occurs. We will seek an optimal position only for the interior nodes; boundary nodes will remain fixed at their original position. To our knowledge this is the first work addressing this issue for the LDG.
The optimization procedure is based on the general framework of a moving mesh for finite element solution of variational problems presented by Tourigny and Hülsemann in [36]. To that end, we first reformulate the LDG discretization as a constrained minimization problem, having as objective function a mesh dependent energy functional. Then, for each grid point, a local variational problem is created in such a way that its energy coincides with the restriction of the global energy functional in that local domain. Instead of using a local Dirichlet problem, as proposed in [36] for the continuous finite element method, an LDG like discretization of a local problem with zero flux boundary conditions is considered. However, in order to guarantee its well posedness a penalization at the boundary of the local domain is imposed. Since the optimal node position is obtained through a local minimization problem no global auxiliary PDE associated to the mesh geometry and monitor function are needed as required by previous work on moving finite element.
We will not strive to seek an explicit representation of the gradient. Instead, a well known gradient free minimization method will be considered for one dimensional problems; and, for two dimensional problems a suitable optimization technique with finite difference approximated gradients will be used. As a result, the algorithm can be applied without major difficulties to high order and even variable degree approximations which is a natural property of all dG methods. In Section 2, we motivate their approach in the context of the projection operator, , where is the natural finite element space, in the dG setting, of discontinuous polynomials.
In Section 3, we briefly remind the LDG method applied to an elliptic equation; and interpret it as a constraint optimization problem. Local variational problems are derived from the global discrete energy functional and are recasted as LDG discretizations of the restriction of the original partial differential equation to small local domains; with zero flux boundary conditions. We show the existence of their solution is guaranteed by imposing a stabilization term at the boundary of the local domain. The optimization procedure applied to the LDG method is summarized in Algorithm 1 . A series of experiments in one and two dimensional domains are presented in Section 4, where the feasibility of the algorithm is shown after only few steps of the optimization procedure.
Section snippets
Motivation of the procedure
Let us briefly motivate the main idea of the procedure on the projection problem. Let be a triangulation of consisting of -simplexes affine equivalent to a reference cell. For every , we denote by the space of polynomials of degree less or equal than . We consider the finite dimensional space ; where no continuity across element boundaries is imposed; and, for simplicity, we have assumed a uniform degree . The orthogonal projection of the function
The Local Discontinuous Galerkin (LDG) method
We consider the following elliptic problem in the domain of , written as a first order linear system: , with Dirichlet boundary conditions, on . For any triangulation of ; consisting of -simplexes affine equivalent to a reference cell; the LDG seeks an approximation , where and is the space of polynomials of total degree less or equal than ; such that the following equations are satisfied on each cell
Numerical experiments
To assess the performance of the -adaptive technique, we have carried out some numerical experiments on a wide variety of one and two dimensional problems with high order approximations. The implementation of the procedure was done in Matlab; and we have extensively used two particular routines of the MatLab Optimization Toolbox: (a) fminbnd which implements the golden ratio method and (b) fmincon for the interior point method. The LDG parameter is in all the tests, however similar
Concluding remarks
An -adaptive procedure based on Tourigny and Hülsemann’s moving mesh algorithm, [36], has been adapted for the LDG spatial discretization applied to a model elliptic problem. To avoid a global minimization problem, which can be difficult to solve, or the use of an auxiliary PDE associated with the mesh geometry and a monitor function as done in the traditional moving finite element approach, the procedure obtains a suitable mesh by (i) considering the LDG formulation as a constraint
References (38)
- et al.
An error indicator monitor function for an r-adaptive finite-element method
J. Comput. Phys.
(2001) - et al.
A theorem for optimum finite-element idealizations
Internat. J. Solids Struct.
(1973) - et al.
The Runge-Kutta discontinuous Galerkin method for conservation laws V: Multidimensional systems
J. Comput. Phys.
(1998) - et al.
An optimal triangulation for second order elliptic problems
Comput. Methods Appl. Mech. Engrg.
(1985) Adaptive grid generation from harmonic maps on Riemannian manifolds
J. Comput. Phys.
(1991)Optimization of finite element grids by direct energy search
Appl. Math. Model.
(1976)Numerical experiments in finite element grids optimization by direct energy search
Appl. Math. Model.
(1977)- et al.
An adaptive moving mesh method with static rezoning for partial differential equations
Comput. Math. Appl.
(2003) - et al.
A two-dimensional moving finite element method with local refinement based on a posteriori error estimates
Appl. Numer. Math.
(2003) - et al.
Moving mesh methods in multiple dimensions based on harmonic maps
J. Comput. Phys.
(2001)