def METHOD_NAME(self, plot=True):
    pde = X2Y2Data()
    mesh = pde.init_mesh(n=0, meshtype='tri')
    space = RaviartThomasFiniteElementSpace2d(mesh, p=0)
    qf = mesh.integrator(3, 'edge')
    bcs, ws = qf.get_quadrature_points_and_weights()
    index = mesh.ds.boundary_edge_index()
    ps = mesh.bc_to_point(bcs, etype='edge', index=index)
    phi = space.METHOD_NAME(bcs, index=index) 
    if plot:
        box = [-0.5, 1.5, -0.5, 1.5]
        fig = plt.figure()
        axes = fig.gca()
        mesh.add_plot(axes, box=box)
        mesh.find_node(axes, showindex=True)
        mesh.find_edge(axes, showindex=True)
        mesh.find_cell(axes, showindex=True)
        node = ps.reshape(-1, 2)
        uv = phi.reshape(-1, 2)
        axes.quiver(node[:, 0], node[:, 1], uv[:, 0], uv[:, 1])
        plt.show()