The paper presents a three-phase algorithm to compute the set of nondominated points for the binary version of the uncapacitated facility location problem with two objectives. The first phase constructs a paving in objective space which is a collection of boxes that covers all nondominated points. The paving procedure is a branch and bound algorithm where boxes are filtered using three dominance-based pruning tests. The second phase is a refinement procedure which attempts to either shrink or eliminate boxes from the paving using supported points calculated in the boxes. The third phase is a generation procedure which solves the allocation problem associated to each of the remaining boxes using a labeling algorithm. The points obtained in each box are filtered over all boxes and provide the set of nondominated points of the initial problem. Computational experiments on standard benchmark instances demonstrate speedups of several orders of magnitude over both a generic algorithm based on a ε-constraint method using Gurobi as MIP solver and a specific algorithm proposed in 2003 by Fernandez and Puerto.