# Interpolate a function given as matrix on a finite element basis

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?

thanks
Antoine

edited Jul 10, 2013

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.