Not yet reviewed by Pith; the record is open.

We demonstrate how to compute the low energy spectrum for small ($N\le 50$), but otherwise arbitrary, spin-glass instances using modern Graphics Processing Units or similar heterogeneous architecture. Our algorithm performs an exhaustive (i.e., brute-force) search of all possible configurations to select $S\ll 2^N$ lowest ones together with their corresponding energies. We mainly focus on the Ising model defined on an arbitrary graph. An open-source implementation based on CUDA Fortran and a suitable Python wrapper are provided. As opposed to heuristic approaches, ours is exact and thus can serve as a references point to benchmark other algorithms and hardware, including quantum and digital annealers. Our implementation offers unprecedented speed and efficiency already visible on commodity hardware. At the same time, it can be easily launched on professional, high-end graphics cards virtually at no extra effort. As a practical application, we employ it to demonstrate that the recent Matrix Product State based algorithm-despite its one-dimensional nature-can still accurately approximate the low energy spectrum of fully connected graphs of size $N$ approaching $50$.