Quiver Inconsistency

mesh
1

I’m not sure if this is merely an artifact of the plotting routine in fenics, but there’s an inconsistency I’d like to resolve. Suppose I have a donut shaped domain and I compute the \hat{\mathbf{r}} unit vector. Using DG0 elements, I would expect these to be located at the cell centers. When I plot the function with plot, I see vectors at vertices, and it’s not entirely clear which ones it is choosing to plot. In contrast, if I plot, manually, using matplotlib's quiver, I get what I expect (see figure). In the figure, I am comparing the color gradient arrows with the black and red arrows. Note that the red arrows are offset to make them more visible; they line up exactly with the black arrows. Essentially, my question is, why do the three sets not line up?

donut

This was generated using the following code:

import numpy as np
from fenics import *
from mshr import *

from matplotlib import pyplot as plt
import sys

r0 = 0.7
r1 = 1.0
Nr = 5

mesh=generate_mesh(Circle(Point(0,0),r1)-Circle(Point(0,0),r0), Nr)
V = FunctionSpace(mesh, "DG",0)
Vvec = VectorFunctionSpace(mesh, "DG",0)
xmesh = V.tabulate_dof_coordinates()[:,0]
ymesh = V.tabulate_dof_coordinates()[:,1]

rvec = Expression(("x[0]/sqrt(x[0]*x[0]+x[1]*x[1])",
    "x[1]/sqrt(x[0]*x[0]+x[1]*x[1])"),element = Vvec.ufl_element())
r_project = np.reshape( project(rvec, Vvec).vector()[:],(47,2))
r_direct = np.asarray([[x/np.sqrt(x**2+y**2),y/np.sqrt(x**2+y**2)] for x,y in zip(xmesh,ymesh)])
plot(mesh)
plt.scatter(xmesh,ymesh)
plt.quiver(xmesh, ymesh, r_project[:,0], r_project[:,1],label="Projection")
plt.quiver(xmesh, ymesh, r_direct[:,0], r_direct[:,1],label="Manual",color='r',pivot='tip')
plot(project(rvec, Vvec))
plt.legend()
plt.savefig("donut.png")
plt.show()
2

dolfin's built in plotting visualises the function at the mesh vertices.

You should use your own plotting algorithm, or output your function using XDMF.write_checkpoint which outputs the function and its function space properties for visualisation with paraview, which handles DG{0,1,2} and CG{1,2} function spaces.

Note that XDMF.write will only output the evaluation of the FE function at the mesh vertices.