Linear scaling integral fitting

In density (or integral) fitting methods, the density (or an orbital product) is replaced with a sum of atom-centered ``auxiliary'' functions, which are used to efficiently compute Coulomb interactions. In this work, we present a method for computing localized fit coefficients that scales linearly with system size, and introduces only extremely modest errors. We apply the algorithm to a variety of methods, including the J piece of the Fock matrix.