Hi, I wanted to know if there is any function to obtain solutions corresponding to the mesh. From the tutorial, if I understanding is correct compute_mesh_vertices() gives the solutions at vertices corresponding to the mesh cell order(mesh.array().). If I project the solutions to a DG0 space, solutions are at the centroid which can be obtained using E.vector().get_local() where E is a functionspace of the obtained solution. Are those in the order of mesh.array() ? how can I get the values for solution corresponding to cell markers which can be obtained from mesh.array().

Here is an illustration showing the 1-1 correspondence between the index of a DG-0 space and the cell index in legacy dolfin:

```
import dolfin
import numpy as np
mesh = dolfin.UnitSquareMesh(10, 10)
V = dolfin.FunctionSpace(mesh, "DG", 0)
u = dolfin.Function(V)
u.interpolate(dolfin.Expression("x[0]+2*x[1]", degree=0))
for cell in dolfin.cells(mesh):
midpoint = cell.midpoint()
np.isclose(u.vector().get_local()[cell.index()], midpoint.x() + 2 * midpoint.y())
```

so if the interpolation is done for the gradient of a function, the gradient at each midpoint is to be found ?

As you are using legacy Dolfin (not DOLFINx), you would have to use a projection.

But yes. If you project the gradient into a DG-0 space, the tabulate_dof_coordinates() of the DG-0 space corresponds to the midpoints of the mesh, with the equivalent values in the `Function.vector().get_local()`

. Note that the easiest way of getting each component of the gradient is to split the projected function with a deepcopy.