# Raviart-Thomas elements on reference square

I'd like to learn how the Raviart-Thomas (RT) element works. To that end I'd like to analytically describe how the basis functions look on the reference square. The goal here is not to implement it myself, but rather just to get an intuitive understanding of the element.

I'm largely basing this work off of the triangular elements discussed here, perhaps extending it to quadrilaterals is a mistake in itself.

That said, I can define the basis functions for first RK element RK0:

$$\mathbf{\phi}_i(\mathbf{x}) = \mathbf{a} + \mathbf{b}\mathbf{x} = \begin{pmatrix} a_1 + b_1 x\\a_2 + b_2 y\end{pmatrix}$$ for $i = 1,\dots,4.$

The conditions on $\mathbf{\phi}_i$ are that:

$$\mathbf{\phi}_i(\mathbf{x}_j)\cdot\mathbf{n}_j = \delta_{ij}$$

where $\mathbf{n}_j$ is the unit normal shown below, and $\mathbf{x}_j$ is its coordinate.

This is the reference square $[-1,1]\times[1,1]$, so this leads to a system of equations for each basis function. For $\mathbf{\phi}_1$ this is:

$$\begin{pmatrix} 1 & 0 & 1 & 0\\ 0 & -1 & 0 & 1\\ -1 & 0 & 1 & 0\\ 0 & 1 & 0 & 1\end{pmatrix}\begin{pmatrix}a_1 \\ a_2\\ b_1\\ b_3\end{pmatrix} = \begin{pmatrix} 1\\0\\0\\0\end{pmatrix}$$

which can be solved to give:

$$\mathbf{\phi}_1(\mathbf{x}) = \frac{1}{2}\begin{pmatrix} 1 + x\\ 0\end{pmatrix}$$

The other basis functions can be found similarly.

Assuming this is correct, the next step is to find the basis functions for RK1. This is where I'm getting a little unsure of myself. According to the link above, the space we are interested in is:

$$P_1(K) + \mathbf{x}P_1(K)$$

A basis for $P_1$ would be $\{ 1, x, y\}$

I think this means the RK1 basis functions should take the form:

$$\mathbf{\phi}_i(\mathbf{x}) = \begin{pmatrix} a_1 + b_1 x + c_1 y + d_1 x^2 + e_1 xy\\ a_2 + b_2 x + c_2 y + d_2 xy + e_2 y^2\end{pmatrix}$$

This leaves 10 unknowns for each basis function. If we apply the same conditions as in the RK0 case, namely:

$$\mathbf{\phi}_i(\mathbf{x}_j)\cdot\mathbf{n}_j = \delta_{ij}$$, where $\mathbf{n}_j$ is the unit normal as shown below:

this gives us 8 equations. The other 2 I think can be found from some moments. I'm not really sure how exactly. The link above talks about integrating against a basis for $[P_1]^2$, but I'm having trouble figuring out what that means. Am I on the right track, or have I completely missed something here?

In general, you cannot just transfer the same polynomial basis from tetrahedral to quadrilateral elements.1 In particular, the whole point of quadrilateral elements is to work with tensor products of one-dimensional polynomials, which is not possible for tetrahedral elements.

There are in fact quadrilateral Raviart-Thomas elements, but their definition is different. In two dimensions, the polynomial space for $RT_k$ is given by $$P_{k+1,k} \times P_{k,k+1},$$ where $$P_{k,l} = \left\{\sum_{i=0}^k\sum_{j=0}^l a_{ij} x^iy^j: a_{ij}\in\mathbb{R}\right\}.$$ A typical polynomial for $k=0$ would thus be as you wrote, but for $k=1$ it would be $$\begin{pmatrix} a_1 + b_1 x + c_1 x^2 + d_1 y + e_1 xy + f_2 x^2y\\ a_2 + b_2 y + c_2 y^2 + d_2 x + e_2 xy + f_2 xy^2 \end{pmatrix}.$$ Hence, $\dim RT_1 = 12$, and in general, $\dim RT_k = 2(k+1)(k+2)$. This means you need two additional degrees of freedom, which should be located on the interior of the element. (In general, for $RT_k$ you take $k+1$ normal derivatives on each facet, and the remaining degrees of freedom from the interior.)

To answer your actual question: For Raviart-Thomas elements, you usually take moments rather than point evaluations, i.e., the remaining conditions come from the conditions $$\int_{-1}^1 \int_{-1}^1 \phi_i(x,y) q_j(x,y) dx\,dx=\delta_{ij},$$ where the $\{q_j\}$ are a basis of $P_{k-1,k}\times P_{k,k-1}$ (e.g., $\{1,x,y\}$ for $k=1$). To make it easier to get a full nodal basis, the facet degrees of freedom are usually not taken as point evaluation, but also as moment conditions: $$\int_{e_m} \phi_i(s)^T\nu_{e_m} \, q_{m,j}(s)\,ds,$$ where $e_m$ is one of the four edges, $\nu_{e_m}$ is the corresponding outer normal, and for each $m$, the $q_{m,j}$ form a basis of $P_k(e_m)$ (e.g., $\{1,x\}$ or $\{1,y\}$ for $k=1$ depending on the edge orientation). Together, these degrees of freedom are unisolvent (i.e., the corresponding system of basis functions is always invertible).

You can find a discussion on quadrilateral Raviart-Thomas-Elements in Chapter 2.4.1 of Boffi, Brezzi, Fortin: Mixed Finite Element Methods and Applications, Springer 2013, Arnold, Boffi, Falk: Quadrilateral $H(div)$ finite elements, SINUM 42(5), 2005, pp. 2429-2451, and Chapter 3.2.3 of Ronald Hoppe's lecture notes.

1. As a rule of thumb, a polynomial space of order $k$ on tetrahedral elements contains monomials whose powers sum to $k$, while a space of order $k$ on quadrilateral elements contains monomials whose maximal power is $k$. For example, $x^2y$ would be of order $3$ on tetrahedrals, but only of order $2$ on quadrilaterals.

• Thank you very much for your answer, you obviously put a lot of effort into it. I think that clears up a lot of my misconceptions. – Lukas Bystricky Jul 22 '15 at 6:12
• I recalculated the basis function $\mathbf{\phi}_1$ for $k=0$ using the integral described above and came up with $\frac{1}{4}\langle 1+x, 0\rangle^T$. Assuming this is correct, could you explain where the compact support comes into play? Since $\mathbf{\phi}_1$ is constant in $y$, it will be non-zero over all the elements above and below it. – Lukas Bystricky Jul 22 '15 at 11:32
• Glad you found it helpful; your question is interesting and you spent a lot of effort as well. The compact support comes from the fact that the polynomials are only defined on the reference element -- recall that Raviart-Thomas are H(div)-conforming elements, and thus functions in the global finite element space need not be continuous. – Christian Clason Jul 22 '15 at 11:34
• Actually, this is only true for the basis functions connected to the interior degrees of freedom: The (global) basis functions connected to the edge degrees of freedom have support on (only) the two elements connected by the edge; on every other element, they are set to zero. – Christian Clason Jul 22 '15 at 11:40
• Actually actually: for edge elements only the normal trace has to be continuous, not the polynomial itself, so even that should be taken care of automatically without extending the support. If you need more details about the global Raviart-Thomas space, I'd suggest you expand your question, and I'll try to expand my answer. – Christian Clason Jul 22 '15 at 13:53