residuals. All operators are implemented using the MeshWorker interface.
Like in step-12, the discretization relies on finite element spaces,
-which are polynomial inside the mesh cells, but have no continuity
-between cells. Since such functions have two values on each face, one
-from each side, we define mean value and jump operators as follows:
-let <i>K</i><sub>1</sub> and <i>K</i><sub>2</sub> be the two cells
-sharing a face, and let the traces of functions <i>u<sub>i</sub></i>
-and the outer normal vectors <b>n<sub>i</sub></i> be labeled
+which are polynomial inside the mesh cells $K\in \mathbb T_h$, but
+have no continuity between cells. Since such functions have two values
+on each interior face $F\in \mathbb F_h^i$, one from each side, we
+define mean value and jump operators as follows: let
+<i>K</i><sub>1</sub> and <i>K</i><sub>2</sub> be the two cells sharing
+a face, and let the traces of functions <i>u<sub>i</sub></i> and the
+outer normal vectors <b>n</b><i><sub>i</sub></i> be labeled
accordingly. Then, on the face, we let
@f[
- \{\!\{ u \}\!\} = \frac{u_1 + u_2}2
+ \{\!\{ u \}\!\} = \frac{u_1 + u_2}2
@f]
Note, that if such an expression contains a normal vector, the
averaging operator turns into a jump. The interior penalty method for the problem
@f[
- -\Delta u = f \;\;\;\; u|_{\partial \Omega} = u^D
+ -\Delta u = f \text{ in }\Omega \qquad u = u^D \text{ on } \partial\Omega
@f]
becomes
-@f[
- \sum_{K\in T_h} (\nabla u, \nabla v)_K
+@f{multline*}{
+ \sum_{K\in \mathbb T_h} (\nabla u, \nabla v)_K
\\
+ \sum_{F \in F_h^i} \biggl\{4\sigma_F (\{\!\{ u \mathbf n\}\!\}, \{\!\{ v \mathbf n \}\!\})_F
- 2 (\{\!\{ \nabla u \}\!\},\{\!\{ v\mathbf n \}\!\})_F
- 2 (\{\!\{ \nabla v \}\!\},\{\!\{ u\mathbf n \}\!\})_F
\biggr\}
\\
- = ???
+ + \sum_{F \in F_h^b} \biggl\{2\sigma_F (u, v)_F
+ - (\partial_n u,v)_F
+ - (\partial_n v,u)_F
+ \biggr\}
+ \\
+ = (f, v)_\Omega + \sum_{F \in F_h^b} \biggl\{
+ 2\sigma_F (u^D, v)_F - (\partial_n v,u^D)_F
+ \biggr\}.
+@f}
+
+Here, $\sigma_F$ is the penalty parameter, which is chosen as follows:
+for a face <i>F</i> of a cell <i>K</i>, compute the value
+@f[
+\sigma_{F,K} = p(p+1) \frac{|F|_{d-1}}{|K|_d},
@f]
+where <i>p</i> is the polynomial degree of the finite element
+functions and $|.|$ denotes the Hausdorff measure of the corresponding
+dimension of the object. If the face is at the boundary, choose $\sigma_F = \sigma_{F,K}$.
+For an interior face, take the average of the two values at this face.
+
+In our finite element program, we distinguish three different
+integrals, corresponding to the sums over cells, interior faces and
+boundary faces above. Since the MeshWorker::loop organizes the sums
+for us, we only need to implement the integrals over each mesh
+element. The class MatrixIntegrator below has these three functions
+for the left hand side of the formula, the class RHSIntegrator for the
+right.
+
+As we will see below, even the error estimate is of the same
+structure, since it can be written as
+@f{align*}{
+ \eta^2 &= \eta_K^2 + \eta_F^2 + \eta_B^2
+ \\
+ \eta_K^2 &= \sum_{K\in \mathbb T_h} h^2 \|f + \Delta u_h\|^2
+ \\
+ \eta_F^2 &= \sum_{F \in F_h^i} \biggl\{
+ 4 \sigma_F \| \{\!\{u_h\mathbf n\}\!\} \|^2 + h \|\{\!\{\partial_n u_h\}\!\}\|^2 \biggr\}
+ \\
+ \eta_B^2 &= \sum_{F \in F_h^b} 2\sigma_F \| u_h-u^D \|^2.
+@f}
+
+Thus, the functions for assembling matrices, right hand side and error
+estimates below exhibit that these loops are all generic and can be
+programmed in the same way.
+
+This program is related to step-12, in that it uses MeshWorker and
+discontinuous Galerkin methods. While there, we solved an advection
+problem, here it is a diffusion problem. Here, we also use mutligrid
+preconditioning and a theoretically justified error estimator, see
+Karakashian an Pascal (2003). The multilevel scheme was discussed in
+detail in Kanschat (2004). The adaptive iteration and its convergence
+have been discussed (for triangular meshes) in Hoppe, Kanschat, and
+Warburton (2009).
\ No newline at end of file