How to use the pyamg.ruge_stuben_solver function in pyamg

# # h = 2.5e-3
    # h = 1.e-1
    # # cell_size = 2 * pi / num_boundary_points
    # c = mshr.Circle(dolfin.Point(0., 0., 0.), 1, int(2*pi / h))
    # # cell_size = 2 * bounding_box_radius / res
    # m = mshr.generate_mesh(c, 2.0 / h)
    # coords = m.coordinates()
    # coords = numpy.c_[coords, numpy.zeros(len(coords))]
    # cells = m.cells().copy()
    # mesh = voropy.mesh_tri.MeshTri(coords, cells)
    # # mesh = voropy.mesh_tri.lloyd_smoothing(mesh, 1.0e-4)

    matrix, rhs = pyfvm.discretize_linear(Poisson(), mesh)

    # ml = pyamg.smoothed_aggregation_solver(matrix)
    ml = pyamg.ruge_stuben_solver(matrix)
    u = ml.solve(rhs, tol=1e-10)
    # from scipy.sparse import linalg
    # u = linalg.spsolve(matrix, rhs)

    mesh.write('out.vtu', point_data={'u': u})
    return

ls = SLS(A=A, b=b)
            sets = self.settings
            sets = {k: v for k, v in sets.items() if k.startswith('solver_')}
            sets = {k.split('solver_')[1]: v for k, v in sets.items()}
            ls.settings.update(sets)
            x = SLS.solve(ls)
            del(ls)
            return x

        # PyAMG
        if self.settings['solver_family'] == 'pyamg':
            if importlib.util.find_spec('pyamg'):
                import pyamg
            else:
                raise Exception('PyAMG is not installed.')
            ml = pyamg.ruge_stuben_solver(A)
            x = ml.solve(b=b, tol=1e-10)
            return x

def _solve_cg_mg(lap_sparse, B, tol, return_full_prob=False):
    """
    solves lap_sparse X_i = B_i for each phase i, using the conjugate
    gradient method with a multigrid preconditioner (ruge-stuben from
    pyamg). For each pixel, the label i corresponding to the maximal
    X_i is returned.
    """
    X = []
    ml = ruge_stuben_solver(lap_sparse)
    M = ml.aspreconditioner(cycle='V')
    for i in range(len(B)):
        x0 = cg(lap_sparse, -B[i].todense(), tol=tol, M=M, maxiter=30)[0]
        X.append(x0)
    if not return_full_prob:
        X = np.array(X)
        X = np.argmax(X, axis=0)
    return X

def _solve_cg_mg(lap_sparse, B, tol):
    if not amg_loaded:
        print """the pyamg module (http://code.google.com/p/pyamg/)
        must be installed to use the amg mode"""
        raise ImportError
    X = []
    #lap_sparse = lap_sparse.tocsr()
    ml = ruge_stuben_solver(lap_sparse)
    M = ml.aspreconditioner(cycle='V')
    for i in range(len(B)):
        x0 = cg(lap_sparse, -B[i].todense(), tol=tol, M=M, maxiter=30)[0]
        X.append(x0)
    X = np.array(X)
    X = np.argmax(X, axis=0)
    return X

# Apply Dirichlet conditions.
        verts = numpy.where(mesh.is_boundary_node)[0]
        # Set all Dirichlet rows to 0.
        for i in verts:
            matrix.data[matrix.indptr[i] : matrix.indptr[i + 1]] = 0.0
        # Set the diagonal and RHS.
        d = matrix.diagonal()
        d[mesh.is_boundary_node] = 1.0
        matrix.setdiag(d)

        rhs = numpy.zeros((n, 2))
        rhs[mesh.is_boundary_node] = mesh.node_coords[mesh.is_boundary_node]

        # out = scipy.sparse.linalg.spsolve(matrix, rhs)
        ml = pyamg.ruge_stuben_solver(matrix)
        # Keep an eye on multiple rhs-solves in pyamg,
        # .
        out = numpy.column_stack([
            ml.solve(rhs[:, 0], tol=tol),
            ml.solve(rhs[:, 1], tol=tol),
        ])
        return out[mesh.is_interior_node]

# Apply Dirichlet conditions.
        verts = numpy.where(mesh.is_boundary_node)[0]
        # Set all Dirichlet rows to 0.
        for i in verts:
            matrix.data[matrix.indptr[i] : matrix.indptr[i + 1]] = 0.0
        # Set the diagonal and RHS.
        d = matrix.diagonal()
        d[mesh.is_boundary_node] = 1.0
        matrix.setdiag(d)

        rhs = numpy.zeros((n, 2))
        rhs[mesh.is_boundary_node] = mesh.node_coords[mesh.is_boundary_node]

        # out = scipy.sparse.linalg.spsolve(matrix, rhs)
        ml = pyamg.ruge_stuben_solver(matrix)
        # Keep an eye on multiple rhs-solves in pyamg,
        # .
        out = numpy.column_stack([
            ml.solve(rhs[:, 0], tol=tol),
            ml.solve(rhs[:, 1], tol=tol),
        ])
        return out[mesh.is_interior_node]

if not isspmatrix_csr(A):
                A = A.tocsr()

            t_solve = time()

            if self.iterative_solver_tolerance > 1e-9:
                self.iterative_solver_tolerance = 1e-10

            # AMG METHOD
            amg_func = None
            if self.preconditioner_type=="smoothed_aggregation":
                # THIS IS TYPICALLY FASTER BUT THE TOLERANCE NEED TO BE SMALLER, TYPICALLY 1e-10
                amg_func = smoothed_aggregation_solver
            elif self.preconditioner_type == "ruge_stuben":
                amg_func = ruge_stuben_solver
            elif self.preconditioner_type == "rootnode":
                amg_func = rootnode_solver
            else:
                amg_func = rootnode_solver

            ml = amg_func(A)
            # ml = amg_func(A, smooth=('energy', {'degree':2}), strength='evolution' )
            # ml = amg_func(A, max_levels=3, diagonal_dominance=True)
            # ml = amg_func(A, coarse_solver=spsolve)
            # ml = amg_func(A, coarse_solver='cholesky')

            if self.solver_context_manager is None:
                # M = ml.aspreconditioner(cycle='V')
                M = ml.aspreconditioner()
                if self.reuse_factorisation:
                    self.solver_context_manager = M