Numerical investigations, error measurements, and convergence tests

Finally, we want to put pieces together. As an orientation we can start with the simple Poisson demos provided in demos and make sure the 1D and 2D demos run through successfully.

Task Solve-1:

Now, we want to measure errors. To this end write a function error(u_h, u_exact) that computes the $L^2$-error between the approximate solution u_h and the exact solution u_exact in the $L^2$ norm, i.e. $$ \int_{\Omega} (u_h - u_{\text{exact}})^2 \, dx. $$ The function should return the error as a scalar. To this end use a loop over elements similar to the one in the assemble-functions. On each element use proper integration rules and proceed similar to the procedure in compute_element_.... functions.

To find a setting where you can verify your implementation, use the following "trick", sometimes known as the method of "manufactured solution":

• Choose a (nice) function $u_{\text{exact}}$ (with proper boundary conditions, i.e. $\partial_n u_{\text{exact}}=0$ on $\partial \Omega$).
• Compute the right hand side $f$ by plugging $u_{\text{exact}}$ into the PDE.
• Solve the PDE numerically with the right hand side $f$ which yields an approximate solution $u_h$.

Task Solve-2:

Next, carry out a convergence study. Let $h$ be the minimal mesh size in your 1D/2D mesh. Observe for a fixed type of finite element space the convergence rate of the error in the $L^2$-norm.

Carry out the experiment for different finite element spaces.

Task Solve-3:

Now, implement the computation of the $H^1$ error: $$ \int_{\Omega} \|(\nabla (u_h - u_{\text{exact}}) \|^2+ (u_h - u_{\text{exact}})^2 \, dx $$ and repeat the convergence studies. What do you observe?

Task Solve-4:

Given a (2D or 1D) domain. We now want to solve the following Poisson problem with homogeneous Dirichlet boundary conditions: Find $u \in H^1_0(\Omega)$ such that $$ \int_{\Omega} \nabla u \cdot \nabla v dx = \int_{\Omega} f v dx $$ for all $v \in H^1_0(\Omega)$.

The discrete problem reads: Find $u_h \in V_{h,0} = \{ v|_T \in \mathcal{P}^k(T) \mid T \in \mathcal{T}_h\} \cap C^0_0(\Omega)$, where $C^0_0(\Omega) := \{v \in C^0(\Omega)\mid v|_{\partial \Omega} = 0\}$ such that $$ \int_{\Omega} \nabla u_h \cdot \nabla v_h dx = \int_{\Omega} f v_h dx $$ for all $v_h \in V_{h,0}$.

To realize this, we first set up the linear system for the unconstrained problem:

Find $u_h \in V_{h} = \{ v|_T \in \mathcal{P}^k(T) \mid T \in \mathcal{T}_h \} \cap C^0(\Omega)$ such that $$ \int_{\Omega} \nabla u_h \cdot \nabla v_h dx = \int_{\Omega} f v_h dx $$ for all $v_h \in V_{h}$. This yields a linear system $A \cdot x = b$ for matrix $A$ and vector $b$.

Then, we need to modify the linear system to account for the Dirichlet boundary conditions.

One simple hack is the following:

Set all rows and columns of $A$ corresponding to degrees of freedom on the boundary to zero, except for the diagonal entries which we set to one. We also set the corresponding entries in $b$ to zero. This yields a linear system $\tilde{A} \cdot x = \tilde{b}$.

Finally, we solve the linear system $\tilde{A} \cdot x = \tilde{b}$ and obtain the solution $u_h$.

Task Solve-5:

Test your implementation of the previous task with a convergence study is in the previous tasks.