Over the past twenty years, the field of plasmonics has been revolutionized with the isolation and utilization of two-dimensional materials, particularly graphene. Consequently, there is significant interest in rapid, robust, and highly accurate computational schemes which can incorporate such materials. Standard volumetric approaches can be contemplated, but these require huge computational resources. Here we describe an algorithm which addresses this issue for nonlocal models of the electromagnetic response of graphene. Our methodology not only approximates the graphene layer with a surface current, but also reformulates the governing volumetric equations in terms of surface quantities using Dirichlet-Neumann Operators. We have recently shown how these surface equations can be numerically simulated in an efficient, stable, and accurate fashion using a High-Order Perturbation of Envelopes methodology. We extend these results to the nonlocal model mentioned above, and using an implementation of this algorithm, we study absorbance spectra of TM polarized plane-waves scattered by a periodic grid of graphene ribbons.