Abstract Finding and characterizing gravitational waves from individual supermassive black hole binaries is a central goal of pulsar timing array experiments, which will require analysis methods that can be efficient on our rapidly growing datasets. Here we present a novel approach built on three key elements: (i) precalculating and interpolating expensive matrix operations; (ii) semi-analytically marginalizing over the gravitational-wave phase at the pulsars; (iii) numerically marginalizing over the pulsar distance uncertainties. With these improvements the recent NANOGrav 15 yr dataset can be analyzed in minutes after an O ( 1 h ) setup phase, instead of an analysis taking days–weeks with previous methods. The same setup can be used to efficiently analyze the dataset under any sinusoidal deterministic model. In particular, this will aid testing the binary hypothesis by allowing for efficient analysis of competing models (e.g. incoherent, monopolar, or dipolar sine wave model) and scrambled datasets for false alarm studies. The same setup can be updated in minutes for new realizations of the data, which enables large simulation studies.