Following the FEniCS Tutorial, I'm trying to retrieve both the nodal coordinates of the mesh, and a function values at those points. This is discussed on pages 17-19 (Section 1.6) of the Langtangen book.
The coordinate values come out just fine, and in the correct order, but the nodal values seem to be in some random order, not what is expected.
Here is some sample Python code:
mesh = UnitSquareMesh(nx,ny)
V = FunctionSpace(mesh, 'Lagrange', degree)
topo = Expression('x[0]')
h = interpolate(topo,V)
coor = mesh.coordinates()
h_nodal_values = h.vector()
something might be going bad here?
h_array = h_nodal_values.array()
if mesh.num_vertices() == len(h_array):
for i in range(mesh.num_vertices()):
print 'h(%4.1f,%4.1f) = %4.1f' % (coor[i][0], coor[i][1], h_array[i])
The output values for a 3 by 3 element mesh, with degree 1, are:
h( 0.0, 0.0) = 1.0
h( 0.3, 0.0) = 1.0
h( 0.7, 0.0) = 0.7
h( 1.0, 0.0) = 1.0
h( 0.0, 0.3) = 0.7
h( 0.3, 0.3) = 0.3
h( 0.7, 0.3) = 1.0
h( 1.0, 0.3) = 0.7
h( 0.0, 0.7) = 0.3
h( 0.3, 0.7) = 0.0
h( 0.7, 0.7) = 0.7
h( 1.0, 0.7) = 0.3
h( 0.0, 1.0) = 0.0
h( 0.3, 1.0) = 0.3
h( 0.7, 1.0) = 0.0
h( 1.0, 1.0) = 0.0
the final numbers after the = sign seem out of order...?