Reduced density matrices of Richardson-Gaudin states in the Gaudin algebra basis

Eigenvectors of the reduced Bardeen-Cooper-Schrieffer Hamiltonian have recently been employed as a variational wavefunction ansatz in quantum chemistry. This wavefunction is a mean-field of pairs of electrons (geminals). In this contribution we report optimal expressions for their reduced density matrices in both the original physical basis and the basis of the Richardson-Gaudin pairs. Physical basis expressions were originally reported by Gorohovsky and Bettelheim. In each case, the expressions scale like $\mathcal{O}(N^4)$, with the most expensive step the solution of linear equations. Analytic gradients are also reported in the physical basis. These expressions are an important step towards practical mean-field methods to treat strongly-correlated electrons.