def METHOD_NAME():
    # Start with a simple field exactly represented in the function space over
    # the unit square domain.
    m = UnitSquareMesh(2, 2)
    V = FunctionSpace(m, "CG", 2)
    x, y = SpatialCoordinate(m)
    f = Function(V).interpolate(x * y)
    # We create a 1D mesh immersed 2D from (0, 0) to (1, 1) which we call "line".
    # Note that it only has 1 cell
    cells = np.asarray([[0, 1]])
    vertex_coords = np.asarray([[0.0, 0.0], [1.0, 1.0]])
    plex = mesh.plex_from_cell_list(1, cells, vertex_coords, comm=m.comm)
    line = mesh.Mesh(plex, dim=2)
    x, y = SpatialCoordinate(line)
    V_line = FunctionSpace(line, "CG", 2)
    f_line = Function(V_line).interpolate(x * y)
    assert np.isclose(assemble(f_line * dx), np.sqrt(2) / 3)  # for sanity
    f_line.zero()
    assert np.isclose(assemble(f_line * dx), 0)  # sanity again
    # We want to calculate the line integral of f along it. To do this we
    # create a function space on the line mesh...
    V_line = FunctionSpace(line, "CG", 2)
    # ... and interpolate our function f onto it.
    f_line = interpolate(f, V_line)
    # The integral of f along the line is then a simple form expression which
    # we assemble:
    assemble(f_line * dx)  # this outputs sqrt(2) / 3
    assert np.isclose(assemble(f_line * dx), np.sqrt(2) / 3)