# discretizing surface integral using nodal DG method

I am currently learning nodal DG methods, primarily through the book by Warburton, and am a bit confused on how to handle surface integrals using straight edged elements. On page 187 (and on page 214) of Warburton's book, we have

$$\int_{\partial E_{k}} \mathbf{n} \cdot f_{h} \ell_{i}^{k}(\mathbf{x}) d\mathbf{x}.$$

Recall that $$f_{h} = \sum_{j=1}^{N} f_{j}(x,t)\ell_{j}(\mathbf{x})$$ and $$\ell_{i}(\mathbf{x})$$ is our basis function, and so,

$$\int_{\partial E_{k}} \mathbf{n} \cdot f_{h} \ell_{i}^{k}(\mathbf{x}) d\mathbf{x} = \int_{\partial E_{k}} \mathbf{n} \cdot \left[ \sum_{j=1}^{N} f_{j}(x,t)\ell_{j}(\mathbf{x}) \right]\ell_{i}^{k}(\mathbf{x}) d\mathbf{x}\\ = \sum_{j=1}^{N} \mathbf{n} \cdot f_{j}(x,t) \int_{\partial E_{k}} \ell{j}(\mathbf{x}) \ell_{i}(\mathbf{x}) d \mathbf{x}\\ = \sum_{j=1}^{N} \mathbf{n} \cdot f_{j}(x,t) M^{k}_{ij}.$$ Now, to discretize the surface mass matrix, we must use a quadrature rule on each face (since I am working in 2D it's a line integral along each edge). With that said, since we have an orthonormal basis, we can compute the mass matrix analytically. For volume integrals, Warburton shows that $$M^{k}_{ij} = (\mathcal{V} \mathcal{V}^{T})^{-1} J,$$ where $$J$$ is the metric Jacobian. So intuitively, the same would be for the mass matrix, however, the Vandermonde matrix, $$\mathcal{V}$$ would be constructed so its on the local face points $$(r,s)$$, where $$(r,s) \in [-1,1]$$. Now, looking at their code (https://github.com/tcew/nodal-dg/blob/master/Codes1.1/CFD2D/CurvedEulerRHS2D.m), we see that the surface integral term has no surface mass matrix. There is an inverse mass matrix applied to the numerical flux but it is the one from the time-derivative. I was curious if someone could point out why there is no surface mass matrix in their formulation.

• I assume that $\ell_j$ is a basis function, but it won't hurt if you clarify your notation. Jun 5, 2021 at 14:25
• For what I understand you I think that you want to interpolate your flux using the same basis functions, that's why you get a mass matrix. But, that's not the only way to solve that integral. You could, for example, use numerical integration and evaluate the function in a discrete set of points. Also, you could assume that the function is constant in your element and use the centroid value. Jun 5, 2021 at 14:28
• I edited my post to specify what l_{j}(x) is. I understand you want to numerically integrate the mass matrix at the quadrature points, however, what I don't understand is why it seems to be 1 in this case (specifically line 56 of the code). In their example code, they interpolate the f_{j} and n to the edges using Vf (Vf = Vg*inv(V), where Vg is the Vandermonde matrix evaluated at the quadrature points). However, why would that remove the need to directly evaluate the mass matrix as the edges? Sorry if I am completely missing something. Wouldnt you interpolate f_{j}, n, and M to the edges? Jun 5, 2021 at 14:48
• You don't need a mass matrix for surface terms. That appears only if you interpolate $f$. I haven't taken a look at the code. Jun 5, 2021 at 15:25

In 2D your volume part consists of a double integral, whereas the surface part is a standard line integral. Generally this is described with tensorial notations.

If you consider the Cartesian case, the volume part results in a mass matrix similar to

$$M_{\Omega}\equiv \mathbf{M}_{\text{1}} \otimes \mathbf{M}_{\text{1}}$$,

whereas the surface part results in a mass matrix similar to

$$M_{\partial \Omega}^{\xi} \equiv \mathbf{M}_{\text{1}} \otimes \mathbf{I}_1\\ M_{\partial \Omega}^{\eta} \equiv \mathbf{I}_1 \otimes \mathbf{M}_{\text{1}}$$.

Here $$\mathbf{I}_1$$ is the identity matrix.

If you invert the volume part and apply it on $$\xi$$-faces you will get something like

$$=\left( \mathbf{M}_{\text{1}} \otimes \mathbf{M}_{\text{1}}\right)^{-1}(\mathbf{M}_{\text{1}} \otimes \mathbf{I}_1) \\ = (\mathbf{M}_{\text{1}} ^{-1} \otimes \mathbf{M}_{\text{1}} ^{-1}) (\mathbf{M}_{\text{1}} \otimes \mathbf{I}_1) \\ = (\mathbf{M}_{\text{1}} ^{-1} \mathbf{M}_{\text{1}} ) \otimes (\mathbf{M}_{\text{1}}^{-1} \mathbf{I}_1)\\ = \mathbf{I}_{\text{1}} \otimes \mathbf{M}_{\text{1}}^{-1}$$

or similar on $$\eta$$-faces

$$=\left( \mathbf{M}_{\text{1}} \otimes \mathbf{M}_{\text{1}}\right)^{-1}(\mathbf{I}_{\text{1}} \otimes \mathbf{M}_1) \\ = (\mathbf{M}_{\text{1}} ^{-1} \otimes \mathbf{M}_{\text{1}} ^{-1}) (\mathbf{I}_{\text{1}} \otimes \mathbf{M}_1) \\ = (\mathbf{M}_{\text{1}} ^{-1} \mathbf{I}_{\text{1}} ) \otimes (\mathbf{M}_{\text{1}}^{-1} \mathbf{M}_1)\\ = \mathbf{M}_{\text{1}}^{-1} \otimes \mathbf{I}_{\text{1}}.$$

Summerizing: I think you are missing the fact, that both mass matrices (volume and surface) cancel, resulting in a 1D inverse mass matrix. Note that numerical fluxes on $$\eta$$-faces are projected into the volume in $$\xi$$-direction, whereas fluxes on $$\xi$$-faces are projected into the volume in $$\eta$$-direction.

Regards