Raviart-Thomas elements for mixed u-p formulation

Question

I'm trying to solve a mixed u-p problem given by the following equations using the finite-element method:

$\int_{\Omega} \delta\epsilon_{ij}2Ge_{ij} \, d\Omega - \int_{\Omega} \delta\epsilon_{ij}\delta_{ij}p \, d\Omega = \int_{\Omega} b_i\delta u_i \, d\Omega + \int_{\Gamma} t_i\delta u_i \, d\Gamma$

$\int_{\Omega} \delta p\left(p/K + \epsilon_{kk} \right) \, d\Omega = 0$

with:

$\epsilon_{ij}=\frac{1}{2}\left(u_{i,j}+u_{j,i}\right)$

$e_{ij}=\epsilon_{ij}-\frac{1}{2}\delta_{ij}p$

$p = -K \epsilon_{kk} = -K u_{k,k}$

and $K$ and $G$ constants.

When I try do discretize $u$ by $u_h \in U_h$ and $p$ by $p_h \in P_h$ such that $U_h \subset H(div, \Omega)$ and $P_h \subset L²(\Omega)$, the solution does not present physical meaning. Depending on the polynomial degrees chosen, I get $NaN$ or invalid oscillations for $u$ or $p$. The same happens if I chose $P_h \subset H¹(\Omega)$.

However, when I try to use $U_h \subset H¹(\Omega)$ and $P_h \subset H¹(\Omega)$ with degrees 2 and 1, respectively, the solution has physical meaning.

By physical meaning I mean that $p$ is constant where it was expected to be constant and $u$ goes to the expected direction in a smooth way.

Question:

Shouldn't $H(div,\Omega)$ be a better space for $u$ since $p$ is proportional to $div(u)$?

Answer 1

To me this looks like the incompressible elasticity equations.

It is well-known in the literature that this problem is analogous to the Stokes flow problem where $P_2-P_1$ is indeed a valid and a stable element as you have noticed. The only difference is the term $1/K (p, \delta p)$, which does not affect the stability proofs in any way.

Shouldn't H(div,Ω) be a better space for u since p is proportional to div(u)?

The space for $u$ and $p$ is determined from the properties of the bilinear form. You simply choose a test function $\delta u = u$ to be equal to the solution and find the term with highest-order derivatives.

In this case, the highest-order term is $$\int_\Omega G \epsilon_{ij} \epsilon_{ij} \,\mathrm{d}\Omega$$ which is finite if and only if $u_i \in H^1(\Omega)$. Thus, $H^1(\Omega)$ is the correct space for $u_i$.

Similar arguments imply that $p \in L^2(\Omega)$. However, since $H^1(\Omega) \subset L^2(\Omega)$, you can sometimes get a stable scheme by further restricting your pressure space and choosing a piecewise-continuous finite element. This is the case in the famous $P_2-P_1$ (aka. Taylor-Hood) element.