How to calculate divergence and vorticity from a velocity field using finite elements

Question

I am in the process of writing a finite element solver for the Navier-Stokes equations and am having trouble computing things like the divergence and vorticity correctly. Currently to compute the vorticity,

$\omega = \frac{\partial{v}}{\partial{x}} - \frac{\partial{u}}{\partial{y}}$,

I would multiply the r.h.s by a test function and integrate over the element. This gives element matrices of the form:

$\Sigma_{i}\int\int_{\Omega}v_{i}\frac{d\phi_{i}}{dx}\phi_{j}dA \hspace{5mm}\text{and}\hspace{5mm} \Sigma_{i}\int\int_{\Omega}u_{i}\frac{d\phi_{i}}{dy}\phi_{j}dA$,

after plugging in our expansions for u and v. This however appears to be the wrong way to do this since my results do not look right. For example after computing the u and v velocity fields:

my vorticity looks like:

In particular the vorticity appears very non-smooth. Here is another example showing the u-velocity and vorticity (zoomed in for clarity):

I have a similar problem with calculating the divergence. It seems that my computed velocity field is not smooth enough and thus when computing the vorticity and divergence I get many positive and negative gradients resulting in choppy fields. So my question is, how do you calculate the vorticity (or divergence) of a computed velocity field using finite elements?

Answer 1

If you have a vector of coefficients $U$ for the velocity approximation, then your velocity field is given by $\vec u_h = \sum_j U_j \vec\varphi_j$. The divergence is then clearly $$ \text{div} \;\vec u_h = \sum_j U_j\;\text{div}\,\vec\varphi_j $$ and similarly for the curl.

The thing to realize, though, is that even if you are using continuous finite elements, i.e., where the $\vec\varphi_j$ are continuous functions, the divergence and curl are not continuous functions. The same then holds true for the divergence and curl of the error. Consequently, you must expect these errors to be oscillatory. This may be what you observe.