  Quantum many-body systems in thermal equilibrium can be described by the
imaginary time Green's function formalism. However, the treatment of large
molecular or solid ab inito problems with a fully realistic Hamiltonian in
large basis sets is hampered by the storage of the Green's function and the
precision of the solution of the Dyson equation. We present a Legendre-spectral
algorithm for solving the Dyson equation that addresses both of these issues.
By formulating the algorithm in Legendre coefficient space, our method inherits
the known faster-than-exponential convergence of the Green's function's
Legendre series expansion. In this basis, the fast recursive method for
Legendre polynomial convolution, enables us to develop a Dyson equation solver
with quadratic scaling. We present benchmarks of the algorithm by computing the
dissociation energy of the helium dimer He$_2$ within dressed second-order
perturbation theory. For this system, the application of the Legendre spectral
algorithm allows us to achieve an energy accuracy of $10^{-9} E_h$ with only a
few hundred expansion coefficients.
