Hi Firedrake-people, I'm working with simple 2D ocean model and I'd like to load bathymetric data to a Function, interpolating to the correct coordinates of the nodes. Currently I get the x/y coordinates from specialized x/y function and then just write the data directly to Function.dat.data array. So something like this: P1 = FunctionSpace(mesh, "CG", 1) bath = Function(P1,name='bathymetry') x_func = Function(P1).interpolate(Expression('x[0]')) y_func = Function(P1).interpolate(Expression('x[1]')) def interpolateBath(x,y): # interpolate here return 3.0 bath.dat.data[:] = interpolateBath(x_func.dat.data, y_func.dat.data) Is there a cleaner way of doing this? I have a similar question for the boundaries. I have external data defined on a line which I'd like to interpolate to the nodes of a certain boundary. I could use a Function for this too, if I knew the indices of the boundary nodes. I've hacked something together using mesh.exterior_facets and FunctionSpace.exterior_facet_boundary_node_map, just wondering if there is an easier way? Regs, Tuomas Karna Post-Doctoral Fellow Oregon Health & Science University