This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

Interpolate a function given as matrix on a finite element basis

+4 votes

Hi all,

Let phi be a function given as a matrix N*N and representing the values of the function at the vertices of a unit square mesh. I want to interpolate this function on the finite element basis. The interpolated function is called psi. I tried this:

mesh = RectangleMesh(0, 0, 1, 1, N, N)
V   = FunctionSpace(mesh, "Lagrange", 1)
psi = Function(V)
phi_vector  = phi.reshape((N+1)*(N+1))
psi.vector()[:] = phi_vector

but the resulting function psi is not correct. I guess this is because the basis functions are not indexed like the mesh vertices. Any idea how to do the interpolation then? Is there a way to find the indexing of the basis functions with respect to the mesh vertices?


asked Jul 9, 2013 by Antoine FEniCS Novice (810 points)
edited Jul 10, 2013 by logg

1 Answer

+6 votes
Best answer

This is the evergreen of the FEniCS questions. There are basically these approaches:

  • enforce ordering of DOFs according to underlying mesh entities (here vertices) by

    parameters['reorder_dofs_serial'] = False

    but this will only work in serial.

  • use helper member functions of V.dofmap(), in this case

    psi.vector()[:] = phi_vector[V.dofmap().vertex_to_dof_map(mesh)]

    should do the job if phi_vector is indexed by local vertex indices.

  • interpolate or project subclassed Expression. This applies well to DG0 spaces as Expression.eval interface supplies cell index to users evaluation algorithm. Hence it is not probably useful for your case.

For further reading refer to

answered Jul 9, 2013 by Jan Blechta FEniCS Expert (51,420 points)
selected Jul 10, 2013 by johanhake

Following your advice I have tried this:

psi.vector()[:] = phi_vector[V.dofmap().vertex_to_dof_map(mesh)]

and it works, many thanks.