How is the integral of a projection over an element $T$ computed in practice? (deal.II related)

Question

I'm studying an error estimator for the equation $\nabla\cdot(\beta u) + cu = f$ and it contains the following term $$||f - cU_h - \Pi(f-c U_h) ||_T$$ where :

• $\Pi$ is the local orthogonal $L^2$ projection onto the classical polynomial DG finite element space $\{v \in L^2: v_{|K} \in \mathbb{P}^{r}(K) \}$
• $U_h$ is the finite element solution, which I have on each cell $T$

I'm implementing this in deal.II. Morally what I want to achieve is the computation of that integral over each cell using quadratures, so

$$\sum_{q} \Bigl(f(x_q) - c U_h(x_q) - \Pi(f-cU_h)(x_q) \Bigr) w_q \det(B_T)$$

I recall that $U_h$ is a vector. Luckily, there's already a routine that allows me to compute $U_h(x_q)$, that is named get_function_values()

However, the true issue is the term $\Pi(f-cU_h)(x_q)$, i.e. the projection of that difference, evaluated at the quadrature points.

What is the practical way people follow to compute that kind of expression? I wasn't able to find anything useful on the deal.II documentation (probably I didn't look in the good places).

---

EDIT:

I tried to follow the second suggestion given by prof. Bangerth. By definition of $L^2$ projection: $$(\Pi(f - cU_h), \phi_i) = (f-cU_h, \phi_i)$$ for every test function $\phi_i$ in my dg space. This means that on every cell I can define and solve $$M_{ij} x_j = (f - c U_h, \phi_i):=F_i$$ where $M_{ij}=(\phi_j,\phi_i)$ is the mass matrix corresponding to a generic cell.

The r.h.s of the equation above is computed by quadrature using the values of $U_h$ at quadrature points (given by the get_function_values() routine), while of course $f$ is a known continuous function evaluated as usual.

I was implementing this by using a MeshWorker framework as done in step-74. That means that the linear system above was put inside the cell_worker lambda expression so that, while looping over all the cells, I compute that projection $x$ for every cell and use that to compute the actual error indicator on each cell.

Answer 1

You need to compute the projection $\Pi(f-cU_h)$ as a first step. This is a finite element field, so you can create a global field (in deal.II lingo: Create the corresponding finite element and a DoFHandler object) and compute it by global projection by inverting a mass matrix. Alternatively, because the field is discontinuous, you can also forgo the global field and just compute the projection locally by inverting a local mass matrix. In either case, you end up with a finite element field on the current cell associated with the discontinuous element, which you can then evaluate at quadrature points like any other finite element field.