We consider the problem of computing a cluster of eigenvalues, and its associated eigenspace, of a (possibly unbounded) selfadjoint operator in a Hilbert space. A rational function of the operator is constructed such that the eigenspace of interest is its dominant eigenspace, and a subspace iteration procedure is used to approximate this eigenspace. The computed space is then used to obtain approximations of the eigenvalues of interest. An eigenvalue and eigenspace convergence analysis that considers both iteration error and discretization error is provided. A realization of the proposed approach for a model second-order elliptic operator is based on a discontinuous Petrov-Galerkin discretization of the resolvent, and a variety of numerical experiments illustrate its performance.