One common argument against the use of discontinuous Galerkin elements
is the large number of globally coupled degrees of freedom that one
-must solve in an implicit system. For the FE_DGP_Monomial basis, each
-scalar solution component is represented by polynomials of degree p
+must solve in an implicit system. This is because, unlike continuous finite
+elements, in typical discontinuous elements there is one degree of freedom at
+each vertex <i>for each of the adjacent elements</i>, rather than just one,
+and similarly for edges and faces. As another example,
+for the FE_DGP_Monomial basis, each
+scalar solution component is represented by polynomials of degree $p$
which yields $(1/dim!)*\prod_{i=1}^{dim}(k+i)$ degrees of freedom per
element. Typically, all degrees of freedom in an element are coupled
to all of the degrees of freedom in the adjacent elements. The resulting
<h4> Reducing the size of the linear system </h4>
To alleviate the computational cost of solving such large linear systems,
-the hybridizable discontinuous Galerkin (HDG) methodology has recently been
-developed by Cockburn and co-workers
+the hybridizable discontinuous Galerkin (HDG) methodology was introduced
+by Cockburn and co-workers
(N.C. Nguyen and J. Peraire:
<i>Hybridizable discontinuous Galerkin methods for partial differential
equations in continuum mechanics</i>, Journal of Computational Physics,
this goal by formulating the mathematical problem using Dirichlet-to-Neumann
mappings. The partial differential equations are first written as a first
order system, and each field is then discretized via a DG method. At this
-point the single-valued ``trace'' values on the skeleton of the
-mesh, i.e. element edges, are taken to be independent unknown quantities.
-The Dirichlet-to-Neumann map concept then permits the solution procedure:
+point the single-valued "trace" values on the skeleton of the
+mesh, i.e. element faces, are taken to be independent unknown quantities.
+The Dirichlet-to-Neumann map concept then permits the following solution procedure:
<ol>
<li> Use local element interior data to enforce a Neumann condition on the
skeleton of the triangulation. The global problem is then to solve for the
(D - C A^{-1} B) \Lambda &=& G - C A^{-1} F \\
A U &=& F - B \Lambda
@f}
-The steps in the Dirichlet-to-Neumann map concept hence correspond to
+The point is that the presence of $A^{-1}$ is not a problem because $A$ is a
+block diagonal matrix where each block corresponds to one cell and is
+therefore easy enough to invert.
+The coupling to other cells is introduced by the matrices
+<i>B</i> and <i>C</i> over the skeleton variable. The block-diagonality of
+<i>A</i> and the structure in <i>B</i> and <i>C</i> allow us to invert the
+matrix <i>A</i> element by element (the local solution of the Dirichlet
+problem) and subtract $CA^{-1}B$ from $D$. The steps in the Dirichlet-to-Neumann map concept hence correspond to
<ol>
<li> constructing the Schur complement matrix $D-C A^{-1} B$ and right hand side $G - C A^{-1} F$,
<li> solving the Schur complement system for $\Lambda$, and
<li> solving the equation for <i>U</i> using the second equation which uses $\Lambda$.
</ol>
-The important ingredient from the linear algebra point of view is that the
-matrix <i>A</i> is block-diagonal with block size equal to the number of
-degrees of freedom of the interior DG variables which are always only related
-to a single cell. The coupling to other cells is introduced by the matrices
-<i>B</i> and <i>C</i> over the skeleton variable. The block-diagonality of
-<i>A</i> and the structure in <i>B</i> and <i>C</i> allow us to invert the
-matrix <i>A</i> element by element (the local solution of the Dirichlet
-problem) and subtract $CA^{-1}B$ from $D$.
<h4> Solution quality and rates of convergence</h4>
Another criticism of traditional DG methods is that the approximate fluxes
and integrate by parts over every element $K$ to obtain:
@f{eqnarray*}
(\kappa^{-1} \mathbf{q}, \mathbf{v})_K - (u, \nabla\cdot\mathbf{v})_K
- + <\hat{u}, \mathbf{v}\cdot\mathbf{n}>_{\partial K} &=& 0, \\
+ + \left<\hat{u}, \mathbf{v}\cdot\mathbf{n}\right>_{\partial K} &=& 0, \\
- (\mathbf{c} u + \mathbf{q}, \nabla w)_K
- + <(\hat{\mathbf{c} u}+\hat{\mathbf{q}})\cdot\mathbf{n}, w>_{\partial K}
+ + \left<(\hat{\mathbf{c} u}+\hat{\mathbf{q}})\cdot\mathbf{n}, w\right>_{\partial K}
&=& (f,w)_K.
@f}
-The terms decorated with a $\hat{}$ denote the numerical traces (also commonly referred
+The terms decorated with a hat denote the numerical traces (also commonly referred
to as numerical fluxes). They are approximations
to the interior values on the boundary of the element. To ensure conservation,
-these terms must be single-valued on any given element edge $\partial K$.
+these terms must be single-valued on any given element edge $\partial K$ even
+though, with discontinuous shape functions, there may of course be multiple
+values coming from the cells adjacent to an interface.
We eliminate the numerical trace $\hat{\mathbf{q}}$ by using traces of the form:
@f{eqnarray*}
\hat{\mathbf{c} u}+\hat{\mathbf{q}} = \mathbf{c}\hat{u} + \mathbf{q}
@f{eqnarray*}
(\kappa^{-1} \mathbf{q}, \mathbf{v})_{\mathcal{T}}
- (u, \nabla\cdot\mathbf{v})_{\mathcal{T}}
- + <\lambda, \mathbf{v}\cdot\mathbf{n}>_{\partial\mathcal{T}}
+ + \left<\lambda, \mathbf{v}\cdot\mathbf{n}\right>_{\partial\mathcal{T}}
&=&
- - <g_D, \mathbf{v}\cdot\mathbf{n}>_{\partial\Omega_D},
+ - \left<g_D, \mathbf{v}\cdot\mathbf{n}\right>_{\partial\Omega_D},
\\
- (\mathbf{c} u + \mathbf{q}, \nabla w)_{\mathcal{T}}
- + <(\hat{\mathbf{c}u} + \hat{\mathbf{q}})\cdot\mathbf{n}, w>_{\partial \mathcal{T}}
+ + \left<(\hat{\mathbf{c}u} + \hat{\mathbf{q}})\cdot\mathbf{n}, w\right>_{\partial \mathcal{T}}
&=&
(f, w)_{\mathcal{T}},
\\
&=&
<g_N, \mu>_{\partial\Omega_N}
@f}
+for all $(\mathbf{v}, w, \mu)$.
<h3> Problem specific data </h3>
- WorkStream to parallelize local solvers. Workstream is already used in step-32, step-44.
- Reconstructing the trace
- Post-processing the solution for superconvergence
-- DataOutFaces: direct output of the global solution
\ No newline at end of file
+- DataOutFaces: direct output of the global solution
+